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Figure  3.1  - 


PPI  photograph  taken  with  10.7  cm  radar  at  3° 
elevation  angle.  Individual  thermals  are  clearly 
seen  as  "doughnut"  shaped  entities  (from  Hardy  and 
Ottersten,  1969). 


Figure  3.2  - 


Figure  3.3  - 


RHI  photographs  at  radar  wavelengths  of  3.2,  10.7, 
and  71.5  cm,  taken  while  the  sky  was  clear. 
Individual  thermals  appear  as  perturbations  near  the 
top  of  the  mixed  layer  at  approximately  1  km  (from 
Hardy  and  Ottersten,  1969. 

Height  contours  of  reflectivity  of  a  typical 
Benard-like  cell.  The  shaded  region  consists  of 
individual  thermals  and  upward  flow,  while  the 
central  region  is  characterized  by  a  slow  downward 
motion  (from  Hardy  and  Ottersten,  1969). 


Figure  3.4  - 


aligned  in 
taken  using 
1.75°  (from 


(a)  An  example  of  individual  thermals  aligned  in 
"thermal  streets".  PPI  photograph  was  taken  using 
S-band  radar  at  an  elevation  angle  of  1.75°  (from 
Konrad,  1970). 

(b)  Contours  of  radar  reflectivity,  C^,  for  a 
two-dimensional  computation  by  Lewellen  (1980).  The 
plot  is  analogous  to  a  RHI  radar  scan  and  shows 
regions  of  high  near  the  top  of  the  mixed  layer, 
Zj.  The  contours  are  labelled  such  that  there  is  a 
factor  of  two  difference  between  each  contour,  and 
the  regions  of  high  would  appear  as  bands  on  a  PPI 
scan. 


Figure  3.5  - 


(a)  Radar  observation  of  a  breaking  Kel vin-Helmholtz 
instability  (from  Gossard,  et  al.,  1971). 

(b)  Reflectivity  contour  of  the  preceeding  radar 
observation. 


Figure  3.6  -  Soliton  solution  of  Eq.  (3.15). 

Figure  3.7  -  Sketch  of  two-sol iton  solution  of  Eq.  (3_J.6).  In 

this  symmetric  case,  Ki  =  K2,  Pi  *  -p2  >  /3KX.  This 
pattern  moves  in  the  x_direction  with  speed 
(Kf  +  pf). 


Figure  3.8  - 


Oblique  interaction  of  two  shallow  water  waves,  taken 
off  the  coast  of  Oregon.  (Photograph  courtesy  of 
Toedtmeier  and  No  Nukes  Calendar,  1979.) 
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Figure  A. 6 


Figure  A. 7 


Figure  A. 8 


(a)  Vertical  distribution  of  virtual  potential 
temperature  during  Flight  #7.  Taken  from  Ruhnke 
(1979). 

(b)  Vertical  distribution  of  specific  humidity  during 
Flight  #7.  Taken  from  Ruhnke  (1979). 

(a)  Vertical  distribution  of  virtual  potential 
temperature  during  Flight  #8.  Taken  from  Ruhnke 
(1979). 

(b)  Vertical  distribution  of  specific  humidity  during 
Flight  #8.  Taken  from  Ruhnke  (1979). 

Isopleths  of  dimensionless  temperature  (T  -  T0)/zJ 
for  the  case  with  Ri  =  0.1.  The  contour  interval  is 
0.2,  and  positive  contours  (i.e.,  lighter  fluid)  are 
denoted  by  a  dashed  line,  (a)  t  =  1.5,  (b)  x  =  2.8, 

(c)  X  =  4.3,  (d)  X  =  5.8,  (e)  X  =  7.3,  and 
(f)  x  =  11.8. 

Isopleths  of  dimensionless  turbulence  kinetic  energy, 
q^/aU^  at  the  same  times  as  Figure  A.l, 
i.e., (a)  x  =  1.5,  (b)  x  =  2.8,  (c)  x  =  4.3, 

(d)  x  =  5.8,  (e)  x  =  7.3,  and  (f)  x  =  11.8.  Contour 
intervals  are  as  follows:  (a)  0.001,  (b)  0.005, 
(c)  0.01,  (d)  0.01,  (e)  0.007,  and  (f)  0.002. 

Evolution  of  total  kinetic  energies,  EK  and  EQ,  for 
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Ri  =  0.1  (.ass. 
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(c)  x  =  9,  (d)  x  =  14.6,  (e)  x  =  20.6.  Contour 
interval  is  0.2. 

Isopleths  of  dimensionless  turbulence  kinetic  energy, 
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(d)  0.003,  (e)  0.001. 

Evolution  of  large  scale  and  small  scale  eddy  kinetic 
energies  for  the  case  with  Ri  =  0.2. 

Profile  of  gradient  Richardson  number,  RiG,  for  the 
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Figure  A. 9  -  Isopleths  of  dimensionless  temperature  variance, 
e^/AT2,  for  the  case  with  Ri  =  0.1.  (a)  t  =  1.5, 

(b)  t  =  2.8,  (c)  t  =  4.3,  (d)  t  =  8.8.  Contour 
intervals  are  (a)  0.004,  (b)  0.02,  (c)  0.02, 

(d)  0.002.  101 

Figure  A. 10  -  Evolution  of  the  total  volume- integrated  temperature 
variance,  ET,  for  the  case  with  (a)  Ri  *  0.1, 

(b)  Ri  =  0.2.  106 

FigureA.il  -  Schematic  representation  of  the  life  cycle  of  an 
individual  Kel vin-Helmholtz  billow  based  on  the  data 
in  the  earlier  figures.  Time  progresses  from  right 
to  left.  Thick  lines  correspond  to  the  detectable 
clear  air  radar  echo,  which  started  as  a  single  layer 
at  1243  and  finished  as  a  double  layer  at  1258  GMT. 

Schematic  vertical  profiles  of  (Ab/Az)  are  indicated 
before  and  after  the  occurrence  of  Kel vin-Helmholtz 


instability.  (From  Browning  and  Watkins,  1970).  109 

Figure  A. 12  -  Evolution  of  eddy  kinetic  energies  for  Ri  =  0.2  with 

(a)  high  resolution  41  x  61,  (b)  low  resolution 
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Figure  A. 13  -  Evolution  of  eddy  kinetic  energies  for  Ri  =  0.2  with 
different  initial  values  of  q2/AU2,  (a)  10-3, 

(b)  3  x  io-3.  113 

Figure  A. 14  -  Evolution  of  eddy  kinetic  energies  for  Ri  =  0.2  with 
different  initial  values  of  A/6,  (a)  0.2,  (b)  0.4, 

(c)  0.8.  116 

Figure  A. 15  -  Evolution  of  eddy  kinetic  energies  for  Ri  =  0.2  with 

different  domain  lengths,  (a)  106,  (b)  156,  (c)  206.  120 
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I.  Introduction 


Modern,  high-powered  radars  can  receive  signals  from  targets  that  are 
nothing  more  than  humidity  and/or  temperature  fluctuations  in  clear  air.  The 
fluctuations  in  refractivity  induced  by  the  fluctuations  in  humidity  or 
temperature  scatter  the  radar  signal  sufficiently  that  the  back  scatter  can  be 
observed  by  a  sensitive  receiver.  Such  clear  air  radar  targets  have  been 
observed  for  some  time  (c.f.,  Ottersten,  Little,  and  Hardy,  1973).  They 
represent  both  an  opportunity  for  probing  atmospheric  motions  and  possible 
interference  with  sensitive  radar  operations.  An  operational,  sensitive 
receiver  must  be  able  to  discriminate  between  these  atmospheric  targets  and 
other  targets  of  potential  interest.  In  some  cases,  unusually  clear  signals 
have  been  associated  with  specific  physical  mechanisms,  such  as  convective 
plumes  (Konrad,  1970)  and  shear-induced  wave  breaking  (Richter,  1979).  More 
frequently  the  precise  nature  of  the  fluctuations  is  insufficiently  known  to 
predict  the  motion  or  evolution  of  the  target.  This  report  presents  the 
results  of  the  initial  phase  of  a  theoretical  investigation  of  different 
possible  fluid  dynamical  mechanisms  which  can  lead  to  relatively  strong  radar 
signals.  The  goal  of  this  study  is  to  provide  a  better  understanding  of  the 
phenomena  and  an  improved  capability  to  predict  the  motion  and  evolution  of 
these  targets. 

This  report  is  divided  into  five  chapters.  Chapter  II  provides  a  brief 
description  of  the  relationships  between  turbulent  fluctuations  in  humidity 
and  temperature  and  the  potential  strength  of  radar  signals.  This  permits  a 
delineation  of  what  types  of  atmospheric  motions  are  detectable  by  radar. 
Chapter  III  provides  a  description  and  rough  cataloguing  of  a  number  of 
different  physical  phenomena  which  may  produce  these  atmospheric  motions. 
Chapter  IV  presents  an  analysis  of  part  of  the  data  obtained  by  NRL  during 
recent  investigations  which  coordinated  radar  probing  with  in  situ  aircraft 
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measurements.  Appendix  A  provides  a  detailed  examination  of  one  particular 
type  phenomena,  Kel vin-Helmholtz  wave  breaking,  using  A. R. A. P.'s  second-order 
closure  turbulence  model.  This  normalized  computation  allows  the  evolution  of 
the  radar  pattern  for  this  particular  type  of  event  to  be  interpreted  in  terms 
of  atmospheric  variables  such  as  the  increments  in  temperature,  humidity,  and 
velocity  across  a  particular  shear  zone,  and  the  wave  speed  of  a  triggering 
internal  wave. 

The  concluding  remarks  provide  specific  recommendations  for  further  work 
which  is  needed  to  extend  our  knowledge  of  the  relationships  between 
particular  features  of  the  radar  return  and  the  nature  of  the  atmospheric 
motion  which  generated  it. 
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II.  The  Dependence  of  Radar  Reflectivity 
on  Humidity  Fluctuations 


Radar  echoes  may  be  generated  by  the  backscatter  caused  by  fluctuations 
in  the  refractive  index.  The  refractive  index  for  cm  wavelength  radio  waves 
is  related  to  temperature  and  humidity  by  the  equation  (Battan,  1970) 

n  =  1  +  10-6^p-  p  +  4810  (2.1) 

with  T  equal  to  temperature  in  °K,  p  equal  to  pressure  in  mb  and  e  equal  to 
water  vapor  partial  pressure  in  mb.  Echoes  generated  by  backscatter  are 
distributed  rather  than  point  targets.  The  radar  equation  for  distributed 
targets  as  given  by  Richter  (1980)  may  be  written  as 

PR  =  0.0177  PTAe  h  n  R-2  (2.2) 

where  Pr  is  the  average  received  power;  Py,  transmitted  power;  Ae,  effective 
antenna  area;  h,  range  resolution;  n,  backscatter  cross  section  per  unit 
volume  (reflectivity);  and  R  is  the  range  of  the  echo  region. 

The  radar  reflectivity  is  related  directly  to  the  spatial  correlation  of 
the  fluctuations  in  refractive  index, 

n  =  0.38  x_1/3  Cp  ,  (2.3) 

where  A  is  the  radar  wavelength,  Cp  is  the  refractive  index  structure  function 
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which  measures  the  spatial  correlation  of  the  fluctuations  in  refracti vity. 


Cn2 


[n(x)  -  n(x  +  r)]2 
r2/3 


(2.4) 


and  r  is  the  separation  distance  for  the  two  point  correlation.  Turbulent 
fluctuations  on  the  centimeter  scale  length  generally  fall  into  a  part  of  the 
inertial  subrange  of  the  turbulent  kinetic  energy  spectrum.  In  this  range  Cj^ 
is  independent  of  separation  distance.  It  is  determined  by  the  total  variance 
in  refractivity  and  the  dissipation  length  scale, 

Cp  «  n' n'  t^2/J  (2.5) 


When  fluctuations  in  pressure  are  neglected  in  comparison  to  fluctuations 
in  temperature  and  water  vapor  mixing  ratio,  the  following  expression  for  Cp 
may  be  derived  (Lewellen  and  Teske,  1975) 


C 


2 

n(radar) 


7.3  *10"3  A‘2/3 


e(,e(, 


26H 


-2(1+  — — *  )(0.61T„  +  B)  e;H' 


+  B2  +  1  + 


2BH 


(O.6IT00)  +  2(0.61Too)B 


H'H' 


(2.6) 


where  B  =  7,730,  p  is  in  atmospheres,  T  is  temperature  in  °K,  H  is  water  vapor 
mixing  ratio  fluctuation,  e(,  is  the  virtual  potential  temperature  fluctuation, 
and  A  is  the  turbulent  macroscale  length  used  in  A.R.A.P.,  Inc.'s  turbulent 
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transport  model . 


In  the  marine  atmospheric  boundary  the  leading  term  in  Cp  will  be 

Cj-  ~  6  *10'b  A"2/jH'H‘  (2.7) 

Thus,  qualitatively,  Cp  can  be  expected  to  be  proportional  to  the  square  of 
the  total  incremental  change  in  water  vapor  mixing  ratio  induced  by  the 
particular  fluid  motion  and  inversely  proportional  to  the  two-thirds  power  of 
the  thickness  of  the  turbulent  zone  across  which  this  change  occurs 
(A  ~  0.2  x  thickness  of  the  turbulent  zone). 

Detectable  values  of  Cp  depend  on  radar  characteristics  and  range  as 
given  in  Eqs.  (2.2  and  2.3)  but  values  as  low  as  10‘15  may  be  significant. 
This  may  be  induced  by  humidity  increments  of  the  order  of  a  fraction  of  a  gm 
per  kg  across  a  turbulent  layer  of  the  order  of  a  few  meters  thick. 

The  task  of  understanding  clear  air  radar  returns  involves  understanding 
particular  fluid  dynamic  features  which  may  be  associated  with  relatively 
strong  changes  in  water  vapor  content  across  relatively  thin  layers  of 
turbulence.  A  catalogue  of  such  features  is  discussed  in  Section  III.  For 
those  features  which  may  be  approximated  as  two-dimensional  in  some  sense,  the 
turbulent  transport  code  developed  at  A.R.A.P.,  Inc.,  can  be  used  to 
quantitati vely  follow  the  evolution  of  the  fluid  dynamic  regions  of 
significant  Cp. 
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III.  A  Catalogue  of  Physical  Phenomena  Leading  to  the 
Necessary  Humidity  Fluctuations  Required  for  Radar  Returns 


Fluid  motions  associated  with  a  number  of  physical  phenomena  can  produce 
clear  air  radar  returns,  and  one  can  catalogue  these  phenomena  using  various 
criteria.  We  initially  choose  two  broad  categories  for  cataloguing,  depending 
on  whether  the  motion  is  associated  with  convective  cells  or  with  internal 
gravity  wave  phenomena. 

3.1  Convective  Motions 

Convective  motion  can  occur  in  situations  where  unstable  density 
gradients  exist.  This  can  occur  either  over  large  horizontal  extents,  such  as 
in  the  atmospheric  boundary  layer  over  land  on  a  sunny  day,  or  in  localized 
regions  of  unstable  gradients  that  may  exist  within  a  stable  background.  An 
example  of  this  latter  situation  is  a  heated  plume  (e.g.,  from  a  power  plant) 
in  the  evening  or  early  morning  atmospheric  boundary  layer,  which  is 
statically  stable.  Theoretical  study  of  such  motions  dates  back  to  the  work 
oi‘  Morton,  Taylor,  and  Turner  (1956),  who  studied  the  motion  of  buoyant 
parcels  using  mixing  length  theory  to  model  the  turbulent  entrainment  of 
ambient  fluid  into  the  parcel .  Although  their  model  is  very  simple,  the  basic 
concepts  used  by  Morton,  et  al.,  are  correct  and  we  now  outline  their  results. 

A  buoyant  parcel  generated  in  an  unstratified  environment  rises  due  to 
gravitational  forces.  Ambient  fluid  is  entrained  into  the  parcel,  causing  a 
dilution  of  the  density  difference  between  parcel  and  environment, 
[pp(z)  -  pamb(z)]>  where  Pp  is  the  parcel  density  and  Pamb  is  the  ambient 
density.  The  dilution  of  the  density  difference  decreases  the  acceleration  of 
the  parcel.  Morton,  et  al.,  showed  that  in  an  unstratified  environment,  the 
vertical  velocity,  w,  parcel  size,  r,  and  density  difference  behave  as 
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w 


~  t"1/2 

r  ~  t1/2 

(pp  "  pamb)  ~  t_3/2  , 

where  t  is  a  measure  of  elapsed  time  from  the  onset  of  the  parcel  motion.  The 
parcel  continues  upward  indefinitely,  although  at  some  point  it  will 
effectively  merge  with  the  environment. 

In  the  case  of  a  stably  stratified  environment,  the  parcel  motion  is 
affected  by  the  ambient  density  gradients.  A  case  of  interest  in  this  report 
is  a  “two-layer"  fluid  which  is  modeled  by  density  pu  overlaying  density  p^, 
where  the  ambient  condition  is  stable,  p^  >  pu,  and  we  define  a  density  jump, 
4p  s  Pi  ■  pu.  A  parcel  released  in  the  lower  region  rises  in  a  manner  similar 
to  that  described  previously  until  it  approaches  the  height  of  the  density 
jump.  For  the  case,  |pu  -  p A | / Ap  >>  1,  the  parcel  is  relatively  unaffected  by 
the  density  jump  and  continues  upward.  The  more  relevant  situation  to  the 
atmospheric  boundary  layer  is  |pu  -  pA|/Ap  <  1.  In  this  case,  the  density 
jump  provides  a  strong  deceleration  of  the  parcel.  Typically,  the  parcel 
slightly  overshoots  the  density  jump  height  and  then  falls  back  and  eventually 
collapses.  The  above  situations  were  studied  in  the  laboratory  by  Scorer 
(1957)  and  the  theory  of  Morton,  et  al.,  agrees  reasonably  well  with 
observations.  It  is  more  difficult  to  apply  the  simple  theories  to  the 
atmospheric  case,  where  the  ambient  turbulence  levels  can  modify  the  parcel 
motion.  In  particular,  the  entrainment  process  can  be  greatly  enhanced  when 
the  environment  is  turbulent,  so  that  the  parcel  decelerates  much  more 
rapidly.  Recent  laboratory  observations  attempting  to  simulate  unstable 
atmospheric  conditions  can  be  found  in  Deardorff,  et  al.  (1980). 
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Thermals  in  the  atmosphere  have  been  extensively  studied  using  techniques 
such  as  glider  planes  (Woodward,  1956),  instrunented  aircraft  (Lenschow  and 
Stephens,  1980)  and  radar  observations  (Hardy  and  Ottersten,  1969;  Konrad, 
1970).  The  following  picture  is  obtained  from  these  studies.  An  individual 
thermal  (referred  to  as  a  "Type  I"  cell  by  Hardy  and  Ottersten)  ranges  in 
horizontal  dimensions  from  about  50  to  100  meters  at  onset  to  about  1000  to 
2000  meters  near  the  top  of  the  mixed  layer,  Zj.  The  parcels  are  warmer  than 
the  environment  up  to  a  level  near  Zj ,  after  which  they  are  cooler  than  the 
surroundings  and  are  negatively  buoyant.  The  parcels  continue  past  their 
neutrally  buoyant  level  due  to  their  upward  momentum.  They  are  likely  to 
possess  a  humidity  excess  throughout  their  ascent.  Turbulent  mixing  and  the 
resulting  refractivity  fluctuations,  Cp,  are  most  pronounced  at  the  outer 
edges,  or  "skin",  of  the  parcel.  This  results  in  clear  air  radar  returns 

which  appear  as  doughnut-shaped  echoes  on  a  "plan  position  indicator"  (PPI) 
scan  (see  Figure  3.1).  These  patterns  become  larger  as  the  height  of  the 
scanning  plane  is  increased.  The  intensity  of  the  radar  returns  depends  on 
the  magnitude  of  C^,  which  is  proportional  to  the  temperature  and  hijnidity 
differences  between  the  thermal  and  the  surroundings,  as  discussed  in 
Chapter  II.  A  "range  height  indicator"  (RHI)  scan  typically  shows  thermals  as 
perturbations  at  the  top  of  the  mixed  layer.  The  mean  hunidity  and 

temperature  gradients  at  the  top  of  the  mixed  layer  give  rise  to  corresponding 
fluctuations  at  this  level.  These  fluctuations  generally  provide  a  radar 
return  which  appears  as  a  continuous  band  centered  near  Zj  on  a  RHI  scan 
(Figure  3.2).  Examples  of  the  radar  returns  from  individual  thermals  can  be 
found  in  Konrad  (1970),  who  used  X-band,  S-band  and  UHF  band  radars  for  his 
studies.  These  radars  had  respective  frequencies  of  9.3  Ghz,  2.8  Ghz  and 
430  Mhz,  and  Konrad  comments  that  the  returns  were  best  for  the  S-band  radar. 
The  thermals  behave  as  independent  entities,  and  as  such,  drift  with  the  local 
mean  wind  as  long  as  there  is  no  strong  vertical  shear  in  the  wind.  Their 

lifetimes  are  of  order  20  to  40  minutes,  which  is  the  time  between  onset  near 

the  statically  unstable  surface  and  the  time  when  the  thermal  reaches  the  top 
of  the  mixed  layer.  For  a  more  detailed  discussion  of  clear  air  convective 
cells  in  the  unstable  boundary  layer,  see  Richter,  Jensen  and  Noonkester 
(1974). 
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Figure  3.1  -  PPI  photograph  taken  with 
10.7  cm  radar  at  3°  elevation  angle. 
Individual  thermals  are  clearly  seen 
as  "doughnut"  shaped  entities  (from 
Hardy  and  Ottersten,  1969). 


Figure  3.2  -  RHI  photographs  at  radar 
wavelengths  of  3.2,  10.7,  and  71.5  cm, 
taken  while  the  sky  was  clear.  Indi¬ 
vidual  thermals  appear  as  perturbations 
near  the  top  of  the  mixed  layer  at 
approximately  1  km  (from  Hardy  and 
Ottersten,  1969). 


Individual  cells  frequently  become  organized  on  a  larger  scale  in  the 
atmospheric  boundary  layer.  There  are  two  ways  in  which  this  collective 
motion  can  happen.  The  first  takes  the  form  of  Benard  cell  motion  when 
observed  in  PPI  scans  (Hardy  and  Ottersten  refer  to  such  motions  as  "Type  II" 
cells).  These  cells  are  typically  5  to  10  km  wide  (as  shown  in  Figure  3.3). 
and  the  overall  structure  can  last  several  hours.  The  outer  boundary  of  the 
Blnard-like  cells  is  composed  of  thermals  with  individual  lifetimes  of  order 
30  minutes,  and  a  mean  updraft  exists  in  this  region.  The  central  area  of  the 
B6nard  cells  is  characterized  by  a  mean  downflow  and  no  clear  air  radar 
returns. 

The  second  type  of  collective  motion  is  commonly  known  as  "thermal 
streets".  In  this  case,  individual  cells  align  parallel  to  the  mean  wind  and 
appear  as  bands  on  PPI  scans.  Thermal  streets  tend  to  appear  when  the  mean 
velocity,  U(z),  is  of  sufficient  magnitude  and  the  wind  direction  is 
sufficiently  constant  to  allow  the  individual  convective  cells  to  become 
aligned.  The  spacing  between  the  bands  ranges  from  approximately  1  to  4  times 
the  mixed  layer  height,  and  che  bands  can  extend  up  to  500  km  in  length.  An 
observational  example  of  thermal  streets  is  shown  in  Figure  3.4a,  and  results 
from  a  two-dimensional  computation  are  presented  in  Figure  3.4b.  It  should  be 
noted  that  the  source  of  the  clear  air  radar  returns  for  the  thermals  need  not 
be  humidity  and  temperature  fluctuations.  Hardy  and  Ottersten  discuss  the 
likelihood  of  insects  acting  as  scatterers.  In  such  cases,  the  insects  are 
caught  in  the  updrafts  and  are  responsible  for  radar  returns. 

3.2  Wave  Phenomena  in  the  Atmosphere: 

Unstable  and  Stable  Disturbances 

We  now  discuss  wave  phenomena  which  can  give  rise  to  clear  air  radar 
returns.  A  distinction  will  be  made  between  (1)  waves  which  are  locally 
generated  and  travel  with  the  local  wind  speed,  and  (2)  waves  which  travel  at 
velocities  different  from  the  local  wind  velocities,  and  propagate  into  the 
region  of  interest  from  afar.  Both  types  of  waves  can  be  understood  using  the 
linearized  Euler  equations,  and  a  brief  theoretical  discussion  of  these  waves 
will  be  presented. 
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Figure  3.3  -  Height  contours  of 
reflectivity  of  a  typical 
Benard-like  cell.  The  shaded 
region  consists  of  individual 
thermals  and  upward  flow,  while 
the  central  region  is  charac¬ 
terized  by  a  slow  downward 
motion  (from  Hardy  and  Otter sten, 
1969). 


i 


i 


20 


Figure  3.4(a)  -  An  example  of  indivi¬ 
dual  thermals  aligned  in  "thermal 
streets".  PPI  photograph  was  taken 
using  S-band  radar  at  an  elevation 
angle  of  1.75°  (from  Konrad,  1970). 


Figure  3.4(b)  -  Contours  of  radar  reflectivity, 
Cp ,  for  a  two-dimensional  computation  by 
Lewellen,  Teske  and  Sheng  (1980).  The  plot  is 
analogous  to  a  RHI  radar  scan  and  shows  regions 
of  high  Cfi  near  the  top  of  the  mixed  layer,  Zj. 
The  contours  are  labelled  such  that  there  is  a 
factor  of  two  difference  between  each  contour, 
and  the  regions  of  high  would  appear  as  bands 
on  a  PPI  scan. 
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The  simplest  configuration  for  studying  atmospheric  waves  is  one  in  which 
the  atmosphere  has  no  regions  of  static  instability,  beginning  at  the  sea 
surface  and  extending  clear  up  to  the  tropopause.  In  this  configuration  all 
small  atmospheric  disturbances  will  be  waves  of  some  type.  To  simplify  the 

presentation,  let  us  assume  temporarily  that  the  atmosphere  is  dry  and 

incompressible.  Also,  we  will  assume  that  the  scale  of  the  motion  is  small 
enough  that  we  may  neglect  the  sphericity  and  rotation  of  the  earth.  The 
motion  of  small  amplitude  disturbances  in  this  configuration  was  analyzed  in 
the  important  work  by  Miles  (1961)  and  Howard  (1961)  on  the  stability  of  a 

sheared  flow  of  a  heterogeneous  fluid  in  the  presence  of  a  constant 

gravitational  field.  Here  we  will  simply  summarize  their  main  results. 

Let  p0(z)  and  U(z)  represent  the  undisturbed  background  distributions  of 
density  and  horizontal  velocity.  Infinitesimal  disturbances  of  the  form 
F(z)  expli k[x-a(k)t]}  are  governed  by  the  linear  differential  equation 

£  [oo(U-a)2^]  +  po  [g*  -  k2(U-a)2jF  =  0,  (3.1) 

where  g  represents  the  gravitational  force  and 


1  dpn 

e(z)  =  - - —  ,  N2  =  gts 

Po  dz 


(3.2) 


N  is  Brunt-Vaisala,  or  buoyancy  frequency,  and  is  the  frequency  (in  rad/sec) 
at  which  a  vertically  displaced  parcel  oscillates.  Equation  (3.1)  is 
sometimes  called  the  Tayl or-Goldstein  equation  (Taylor,  1931;  Goldstein, 
1931).  When  coupled  with  appropriate  boundary  conditions  at  the  sea  (or  land) 
surface  and  at  the  tropopause  [e.g.,  one  choice  is  F(0)  *  F(H)  =  0,  since  the 
vertical  velocity  is  proportional  to  F]  it  may  be  considered  a  linear 
eigenvalue  problem  for  a(k),  the  eigenvalue.  Miles  and  Howard  proved  that: 

(i)  the  flow  is  not  unstable  if 
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i.e.,  the  local  Richardson  number  must  be  less  than  0.25  somewhere 
in  the  flowfield  for  instability. 

(ii)  If  there  are  any  unstable  modes,  with  complex  phase  speed 

a  =  ar  +  i  a-j  ,  a^  >  0  (3.4) 

then 

a  =  min[U(z)]  <  ar  <  max[U(z)]  =  n  ,  (3.5) 

z  z 

and 

j^ar-|(a  +  e)j  +  af  <  |'(o-e)J  .  (3.6) 

The  last  result  is  known  as  the  semi -circle  theorem. 

3.2.1  Unstable  atmospheric  disturbances.  It  follows  from  Eq.  (3.5)  that 
associated  with  each  unstable  wave  mode  is  at  least  one  so-called  critical 
layer,  where  U(z*)  =  ar;  i.e.,  where  the  wave  travels  with  the  wind.  The 
unstable  modes  are  collectively  referred  to  as  stratified  shear- flow 

instabilities  when  the  density  gradient  provides  a  stable  stratification, 
i.e.,  e  >  0.  The  most  extensively  studied  mode  of  instability  is  that 

associated  with  an  inflection  point  in  the  velocity  profile.  This  type  is 

probably  also  the  most  prevalent  stratified  shear-flow  instability  in  the 
atmosphere,  and  is  frequently  referred  to  as  Kel vin-Helmholtz  instability, 
although  the  latter  is  strictly  confined  to  discontinuous  velocity  and 
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temperature  profiles.  In  this  report,  we  shall  use  the  phrase 
Kelv in-Helmholtz  (K-H)  instability  to  denote  the  shear-induced  instability 
associated  with  a  confined  layer  of  velocity  and  temperature  gradient.  A 
detailed  study  of  the  K-H  instability  is  given  in  Appendix  A.  Equation  (3.3) 
states  that  for  unstable  waves  to  exist,  Ri  5  g ti/ ( dU/dz ) 2  <  0.25  is  a 
necessary,  but  not  sufficient,  condition.  This  condition,  Ri  <  0.25,  has  the 
physical  meaning  that  the  destabilizing  effect  of  the  mean  shear  must  be 
greater  than  the  stabilizing  effect  of  the  density  stratification.  In 
general,  shear-induced  waves  are  found  in  atmospheric  regions  where  Ri  <  0.25. 
Kel vin-Helmhol tz  instabilities  commonly  occur  at  high  altitudes  in  the 
atmosphere  (5  to  7  km)  near  the  jet  stream  where  velocity  shears  are  large, 
and  in  the  boundary  layer  (~  200  to  1000  meters)  where,  again,  one  finds  large 
shears.  As  this  report  is  primarily  concerned  with  clear  air  radar  returns 
originating  in  the  atmospheric  boundary  layer,  only  a  brief  discussion  of  K-H 
waves  in  the  upper  atmosphere  will  be  given. 

Reed  and  Hardy  (1972)  and  Hooke  and  Hardy  (1975)  combined  measurements 
from  a  sensitive  3  Ghz  radar  with  surface  microbarographs  and  conventional 
synoptic  data  to  document  a  stratified  shear-flow  instability  which  occurred 
over  the  eastern  seaboard.  The  waves  possessed  a  wavelength  and  frequency  of 
approximately  15  km  and  300  seconds,  respectively,  and  a  horizontal  phase 
velocity  of  about  45  m/sec.  The  instability  was  associated  with  the  high 
shear  region  near  the  lower  portion  of  the  jet  stream  (z  ~  8  km),  and  a 
critical  layer  was  found  at  z  =  9  km.  Radar  returns  on  a  RHI  scan  revealed  a 
number  of  echo  strata  throughout  heights  of  8  to  10  km.  The  returns  appeared 
as  undulating  strata  and  returns  near  the  critical  level  showed  the  "cat's 
eye"  pattern  typical  of  K-H  waves,  as  shown  in  Figure  3.5.  The  waves  are 
vertically  trapped  with  amplitudes  decaying  away  from  the  critical  level; 
typically  velocity  fluctuations  are  negligible  1  km  above  or  below  the 
critical  layer.  Two  more  recent  studies  of  unstable  waves  in  the  upper 
atmosphere  were  carried  out  by  VanZandt,  Green,  Clark,  and  Grant  (1979),  and 
Klostermeyer  and  Ruster  (1980),  using  doppler  radar.  The  dynamics  of  the 
waves  observed  by  them  were  similar  to  that  just  described,  but  both  noted  the 
feature  of  "power  bursts".  These  bursts  were  due  to  enhanced  radar  returns 
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Figure  3.5(b)  -  Reflectivity  contour  of  the  preceeding 
radar  observation. 
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from  localized  regions  of  wave  overturning.  That  is,  the  instability  grew  to 
the  point  where  heavy  fluid  was  advected  over  light  fluid;  this  statically 
unstable  situation  then  generated  high  levels  of  turbulence  and  Cp  locally. 
This  development  will  be  demonstrated  in  the  model  calculations  of  Appendix  A. 

Next  we  discuss  wave  phenomena  in  the  atmospheric  boundary  layer  which 
can  give  rise  to  clear  air  radar  returns.  Pioneering  work  in  probing  the 
lower  atmosphere  with  high  resolution  radar  was  done  by  Richter  and 
co-workers,  who  developed  a  high-sensitivity,  high  resolution  radar  system  at 
the  Naval  Electronics  Laboratory  Center  in  San  Diego.  The  radar  operated  at 
2.9  Ghz  and  radar  measurements  were  made  along  the  San  Diego  coast.  The 
atmospheric  boundary  layer  in  this  region  is  characterized  by  a  stable 
inversion,  typically  above  200  to  300  meters,  and  wind  speeds  of  order  1 
m/sec.  The  principal  results  of  Richter  and  co-workers  studies  are  given  in 
three  papers:  Gossard,  Richter  and  Atlas  (1970);  Atlas,  Metcal i ,  Richter  and 
Gossard  (1970);  and  Gossard,  Jensen  and  Richter  (1971).  These  papers  treat 
various  phenomena  responsible  for  clear  air  radar  returns,  and  we  will  first 
discuss  shear  instabilities.  These  instabilities  typically  occurred  near  the 
base  of  the  inversion  (i.e.,  near  the  top  of  the  mixed  layer),  ~  200  to  400 
meters,  and  temperature- velocity  soundings  showed  this  region  to  be  one  with 
critical,  or  near  critical,  Richardson  numbers.  The  waves  appeared  on  a  RHI 
scan  as  well-defined  rolls  with  peak  to  peak  amplitudes  of  up  to  100  meters. 
Typical  wavelengths  and  periods  were  200  meters  and  250  seconds,  respectively, 
and  phase  velocities  were  of  order  1  m/sec.  Again,  critical  layers  were 
present  where  the  phase  velocity  equaled  the  mean  velocity  (z  ~  200  to  400 
meters).  Superimposed  on  the  rolls  were  small  scale  perturbations  due  to 
overturning  (localized  regions  of  static  instability).  These  perturbations 
represent  the  same  phenomena  which  gave  rise  to  the  "power  bursts"  in  the 
study  of  upper  atmospheric  K-H  waves.  Evidence  of  long  period  internal 
gravity  waves  was  found  in  all  of  the  above  studies,  and  it  was  proposed  that 
these  waves  were  instrumental  in  triggering  the  shear  instabilities;  the 
mechanism  is  as  follows.  The  base  of  an  inversion  (or  top  of  a  mixed  layer) 
typically  has  a  mean  velocity  shear.  Situations  can  exist  where  the 
Richardson  number  is  marginally  stable  (Ri  slightly  above  0.25),  and  the  shear 


of  a  propagating  internal  wave  can  add  to  the  background  shear,  driving  the 
Richardson  nunber  below  0.25.  We  return  to  this  scenario  in  Section  3.2.2. 

A  number  of  other  observations  of  waves  in  stratified  shear  flows  can  be 
found  in  the  literature.  Noonkester  (1973)  compiled  data  on  a  total  of  44 
wave  events  which  occurred  over  the  Southern  California  coastal  region. 
Hooke,  Hall  and  Gossard  (1973)  observed  waves  over  Colorado  using  acoustic 
sounders  and  an  instrumented  tower.  Generally  in  boundary  layer  studies, 
wavelengths  are  on  the  order  of  several  hundred  meters,  which  is  about  6  to  10 
times  the  thickness  of  a  typical  shear  layer,  and  phase  velocities  are  a  few 
meters  per  second,  corresponding  to  a  typical  mean  wind  velocity.  A  fairly 
complete  list  of  observations  of  shear  instabilities  can  be  found  in  a  review 
article  by  Ottersten,  Hardy  and  Little  (1973). 

We  also  note  that  an  interaction  between  the  convective  cells  considered 
in  Section  3.1  and  the  stable  layer  capping  the  boundary  layer  can  also  give 
rise  to  K-H  instabilities.  When  there  is  significant  wind  shear  in  this 
stable  layer  as  is  often  the  case,  the  individual  convective  cell  protruding 
up  into  this  layer  forms  a  hummock  with  increased  shear  across  the  top  of  it. 
K-H  wave-breaking  has  also  been  observed  at  the  top  of  these  hummocks 
(Browning,  et  al .  1973).  Due  to  the  vertical  wind  shear  the  resulting  regions 
of  enhanced  Cp  will  move  at  a  velocity  different  from  the  wind  velocity  in  the 
mixed  layer. 

3.2.2  Stable  wave  modes.  Equation  (3.5)  states  that  in  the  absence  of 
external  forcing,  any  wave  whose  speed  differs  from  the  wi nd  everywhere  (i.e., 
for  all  z)  cannot  be  unstable.  This  allows  atmospheric  disturbances  to  be 
divided  into:  (i)  plumes  of  various  sorts,  all  of  which  are  carried 
(horizontally)  by  the  wind;  (ii)  unstable  traveling  waves,  each  of  which 
travels  with  the  wind  at  its  own  critical  level;  and  (iii)  neutrally  stable 
traveling  waves,  which  do  not  travel  with  the  wind  at  any  level.  Categories 
(i)  and  (ii)  have  been  discussed  in  our  treatment  of  thermals  and 
shear-induced  waves,  respectively.  The  stable  atmospheric  disturbances  of 
category  (iii)  will  now  be  considered. 
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The  derivation  of  Eq.  (3.1)  was  based  on  the  assumptions  that  the 
atmosphere  is  incompressible  and  dry.  We  may  relax  these  assumptions  if  we 
require  that: 

(i)  the  motion  is  isentropic; 

(ii)  the  specific  humidity,  q,  of  each  fluid  particle  is  conserved  during 
the  motion  (i.e.,  no  condensation  or  evaporation  occurs); 

(iii)  the  moist  air  is  a  perfect  gas,  with  gas  constant 


R  =  Ra  [1  +  0.61  q0(z)] 


(3.7) 


where  Rd  is  the  gas  constant  for  dry  air,  and  q0(z)  is  the 
background  distribution  of  specific  humidity  (e.g.,  see  Holmboe, 
Forsythe,  and  Gustin,  1943). 

It  is  convenient  to  define  the  virtual  temperature 


T$(z)  »  [1  +  0.61  q0(z)]T0(z)  , 


(3.8) 


where  T0(z)  is  the  background  temperature.  Then  an  analysis  similar  to  those 
of  Miles  and  Howard  leads  again  to  Eq.  (3.1),  provided  we  replace  t{z)  in 
Eq.  (3.2)  with 


1  1 

'£$  +  _£_| 

/I  +  0.61  q0\" 

T* 

'0  i 

_  dz  Cp<3 

^1  +  0.90  q0J 

(3.9) 


where  Cpa  is  the  specific  heat  of  dry  air  at  constant  pressure.  As  before, 
static  stability  requires  e  >  0,  while  dynamic  stability  is  assured  if 
Eq.  (3.3)  holds.  Moreover,  it  follows  again  from  Eq.  (3.5)  that  only 
neutrally  stable  waves  can  travel  at  speeds  different  from  the  wind. 
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For  convenience,  we  have  represented  the  tropopause  (at  z  =  H)  by  a  rigid 
lid  [i.e.,  F(H)  =  0],  However,  no  qualitative  changes  occur  if  the  tropopause 
is  represented  by  a  free  surface.  In  particular,  it  remains  true  that  in  the 
absence  of  external  forcing,  only  neutrally  stable  gravity  waves  can  travel  at 
speeds  different  from  the  wind  at  all  levels.  This  is  a  fundamental  feature 
of  propagating  gravity  waves. 

As  examples,  let  us  consider  two  special  cases  of  cs(z)  and  tl(z)  which 
have  application  to  the  atmosphere.  The  first  is  the  case  of  constant  tf(z), 
and  we  take  U  =  0  for  simplicity.  The  solutions  to  Eq.  (3.1)  are 


F  +  =  F0  exp(±imz)  ,  (3.10) 

which  represent  freely  propagating  disturbances.  In  the  case  of  rigid 
boundaries  at  z  =  0  and  z  =  H,  a  linear  combination  of  the  two  solutions 
satisfies  the  boundary  conditions.  The  dispersion  relation  for  these 
disturbances  is 


Qti  =  Nz  k^ 

k2  +  mz  k2  +  mz 


(3.11) 


from  which  it  follows  that  | u> |  <  N  and  that  |uj|  ♦  N  as  m/k  ♦  0,  | «*»|  +  0  as 
k/m  +  0  display  this.  One  may  show  that  horizontally  propagating  waves  are 
long  waves  (k  <<  m) ;  see  Mowbray  and  Rarity  (1967). 

A  second  form  of  d  applicable  tc  the  atmosphere  is  one  which  models  an 
atmosphere  often  observed  over  the  ocean,  and  similar  to  the  environment 
studied  by  Richter  and  co-workers  at  San  Diego.  The  region  from  approximately 
200  to  1000  meters  is  very  stably  stratified,  while  the  region  above  1000 
meters  is  less  so.  We  take  the  following  form  for  N(z)  =  [gp(z)J1/z, 
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200  m  <  z  <  1000  m  , 


N(z)  =  Nl 

N(z)  =  Nu  1000  m  <  z  , 


where  Nj_  >  Nu  and  the  subscripts  denote  lower  and  upper,  respectively.  As 
above,  we  assume  U(z)  =0.  A  set  of  eigenfunctions  exist  with  Nu  <  w  <  NL. 
These  solutions  oscillate  in  the  region  200m  <  z  <  1000m,  and  decay 
exponentially  for  z  >  1000m.  Such  waves  are  "ducted"  modes.  A  set  of 
eigenfunctions  also  exists  with  id  <  Nu,  and  these  disturbances  can  freely 
propagate  throughout  the  atmosphere.  Both  types  of  waves,  freely  propagating 
and  ducted,  are  commonly  found  in  the  atmosphere.  The  former  are  low 
frequency  disturbances  [uo  <  N(z)  through  tie  atmosphere],  while  the  latter  are 
high  frequency  modes  which  cannot  propagate  into  regions  where  <d  >  N(z). 

The  preceeding  discussion  excluded  any  mean  velocity.  However,  in  the 
case  of  constant  U(z),  it  follows  from  Eq.  (3.1),  that  the  effect  of  the  mean 
flow  is  merely  to  "doppler-shift"  the  frequency  of  the  wave  relative  to  a 
stationary  observer.  The  problem  becomes  more  difficult  to  solve  when  shear 
exists  (aU/az  i  0),  but  for  cases  where  a  >>  U,  the  effects  of  U(z)  will  not 
be  dramatic.  On  the  other  hand,  if  U(z)  =  a  somewhere  in  the  flow,  we  may 
anticipate  the  unstable  waves  previously  discussed. 

The  existence  of  propagating  internal  gravity  waves  in  the  atmospheric 
boundary  layer  has  been  documented  by  Gossard,  Jensen  and  Richter  (1971).  In 
certain  cases,  multi-layered  structures  appeared  on  RHI  scans  which  could  be 
interpreted  as  marking  the  crests  of  a  propagating  wave.  As  the  clear  air 
returns  marked  regions  of  enhanced  ,  Gossard,  et  al .  (1971),  investigated 
mechanisms  which  might  provide  increased  turbulence  for  a  propagating  internal 
wave.  The  first  mechanism  investigated  was  that  of  wave  overturning,  where 
the  amplitude  of  the  wave  is  large  enough  to  advect  heavy  fluid  over  light, 
thereby  creating  localized  regions  of  static  instability.  Data  did  not 
support  this  mechanism  as  balloon  soundings  showed  no  evidence  for  thermal 
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(convective)  instability  at  the  heights  of  radar  returns.  Gossard,  et  al., 
did  find  that  a  high  resolution  (of  order  of  meters)  inspection  of  the  echo 
bands  revealed  a  sub-structure  remarkably  similar  to  that  seen  in  large  scale 
shear-induced  instabilities.  Gossard,  et  al.,  proposed  that  internal  waves 
propagate  into  a  region  possessing  sub-critical  ambient  shear,  and  that  the 
sum  of  wave  shear  plus  ambient  shear  creates  a  dynamically  unstable  situation 
(Ri  <  0.25).  They  argued  that  the  conditions  for  this  dynamic  instability 
will  occur  before  the  onset  of  a  convective  instability  and  the  dynamic 
instability  are  responsible  for  the  layered  radar  returns.  These  arguments 
were  refined  and  extended  in  a  later  work  by  Gossard,  Richter  and  Jensen 
(1973).  We  believe  that  this  mechanism  is  a  prime  candidate  for  clear  air 
radar  returns  which  propagate  at  velocities  significantly  different  from  the 
local  wind  velocity.  Internal  waves  typically  can  have  phase  velocities  much 
greater  than  the  local  wind  speed  and  arbitrary  directions  with  respect  to  the 
local  wind.  They  can  create  regions  of  critical  Ri  which  propagate  at  their 
phase  speed,  giving  rise  to  propagating  regions  of  dynamic  instability  and 
enhanced  Cj^.  As  long  as  the  characteristic  time  for  the  evolution  of  the  wave 
breaking  generated  in  this  unstable  region  is  smaller  t,.an  the  time  required 
for  the  wave  to  travel  through  the  region,  the  region  of  enhanced  will 
appear  to  propagate  with  the  speed  of  the  wave.  This  evolution  of  the 
breaking  wave  is  studied  in  detail  in  Appendix  A.  Given  the  demonstrated 
ability  of  internal  gravity  waves  to  produce  clear  air  returns,  and  the 
proposed  way  in  which  internal  waves  can  give  rise  to  rapidly  moving  radar 
returns,  we  proceed  to  discuss  the  observations  and  production  mechanisms  for 
internal  gravity  waves  in  the  atmosphere. 

Freely  propagating  gravity  waves  have  been  observed  by  a  number  of 
workers,  and  Jordan  (1972)  studied  waves  generated  by  winds  and  storms  in  the 
Colorado  area  using  a  ground  based  array  of  microbarographs  and  anemometers. 
Certain  wave  phenomena  were  related  to  topography,  such  as  mountain  lee  waves, 
while  other  pressure  oscillations  were  due  to  shear  instabilities  associated 
with  either  the  jet  stream  or  the  shear  due  to  frontal  systems.  In  addition, 
Jordan  observed  internal  gravity  waves  generated  from  summer  thunderstorms 
with  periods  of  order  20  to  30  minutes,  wavelengths  of  order  30  km,  and 
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horizontal  phase  speeds  of  order  30  m/sec.  Jordan  also  observed  pressure 
pulses  on  the  microbarographs  which  were  followed  by  nigh  frequency  (periods 
of  a  few  minutes)  oscillations.  Two  mechanisms  for  wave  generation  were 
offered  by  Jordan.  The  first  proposed  was  that  the  downdraft  cells  of  a 
thunderstorm  impinge  on  the  boundary  layer  inversion  which  perturbs  the  region 
and  excites  internal  waves.  The  second  mechanism  involves  the  updraft  region 
of  the  thunderstorm,  where  fluid  can  be  carried  above  its  equilibrium 
position.  Fluid  parcels  will  "overshoot"  equilibrium  levels  and  perturb  the 
surroundings,  which  will  excite  internal  waves.  A  visual  display  of  internal 
gravity  waves  due  to  a  severe  thunderstorm  can  be  found  in  the  satellite 
pictures  of  Erickson  and  Whitney  (1973).  Bands  of  internal  waves  propagate 
out  from  a  storm  located  in  the  Southern  Great  Plains.  The  wavelength  is 
approximately  5  to  10  km,  and  the  clouds  are  estimated  to  be  at 
mid-tropospheric  levels  (~  4  to  5  km),  so  that  one  has  a  picture  of  internal 
gravity  waves  extending  throughout  the  troposphere  and  propagating  away  from 
the  storm.  Uccelini  (1975)  used  synoptic  data  to  infer  the  existence  of  very 
long  wavelength,  low  frequency  waves  associated  with  a  severe  midwest  storm. 
Wavelengths  and  periods  were  approximately  400  km  and  2  hours,  respectively, 
which  corresponds  to  a  phase  velocity  of  approximately  40  m/sec.  Gossard, 
Richter  and  Atlas  (1970)  observed  internal  gravity  waves  on  the  Gan  Oiego 
coast  using  radar  and  ground  based  microbarographs  and  anemometers.  Wave 
periods  ranged  from  7  to  20  minutes,  and  very  large  amplitudes  were  observed. 
An  echo  layer  at  approximately  200  meters  showed  peak  to  peak  amplitudes  of  up 
to  100  meters.  The  wave  form  was  highly  distorted  due  to  non-linear  effects 
with  steepened  crests  and  flattened  troughs.  No  evidence  of  wave  breaking  was 
found  even  with  this  high  degree  of  non-linearity.  Wavelengths  and  phase 
velocities  ranged  from  3  to  7  km  and  2  to  7  m/sec,  respectively.  These  phase 
speeds  always  exceeded  the  mean  winds,  so  that  no  evidence  was  found  for 
critical  layers.  Gossard,  et  al.,  also  speculated  on  the  wave  source  and 
noted  that  thunderstorms  were  present  some  500  km  from  their  instruments. 
More  recent  studies  have  added  to  the  evidence  for  storm  generated  internal 
waves.  Oliver  (1981)  discusses  satellite  data  taken  over  the  eastern  Pacific 
which  clearly  shows  waves  originating  in  the  vicinity  of  strong  fronts  and 
propagating  toward  the  equator.  Gedzelman  (1981)  analyzed  wave  events  at 
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Palisades,  New  York,  covering  a  one-year  period.  Approximately  one-half  of 
these  events  could  be  traced  to  shear  instabilities.  Most  of  the  others  were 
found  to  propagate  at  velocities  significantly  different  from  the  mean  wind 
aloft  and,  in  some  cases,  could  be  traced  to  thunderstorm  activity. 

To  summarize,  thunderstorms  are  a  common  source  of  freely  propagating 
gravity  waves  in  the  atmosphere.  They  are  by  no  means  the  only  source  of 
these  waves.  In  the  lee  of  large  mountain  ranges,  such  as  the  Rockies  and  the 
Sierras,  large  amplitude  lee  waves  are  common  and  are  well  documented  (e.g., 
Lilly  and  Toutenhoofd,  1969;  Lilly,  1971).  On  a  smaller  scale,  lee  waves 
also  may  be  generated  in  the  lee  of  heat  sources,  such  as  a  city  over  land 
(e.g.,  Wurtele,  1955)  or  an  island  in  the  ocean.  One  mechanism  for  generation 
that  is  relatively  rare  was  identified  by  VanDorn  (1964).  He  noted  that  the 
deformation  of  the  earth's  surface  by  the  Alaska  earthquake  of  1964  created  an 
"atmospheric  tsunami"  that  was  observed  as  far  away  as  northern  California. 

The  existence  of  ducted,  or  trapped,  internal  waves  has  been  well 
documented  by  Metcalf  (1975),  who  used  an  instrumented  aircraft  in  conjunction 
with  radar  and  balloon  soundings  to  study  waves  near  the  San  Diego  coast. 
Radar  returns  showed  that  the  low  level  inversion  oscillated  almost 
continuously,  while  analysis  of  aircraft  temperature  and  velocity  records 
revealed  the  presence  of  internal  waves  with  wavelengths  and  periods  ranging 
from  40  to  700  meters  and  0.50  to  1.50  minutes,  respectively.  Nearly  all  of 
the  observed  waves  were  trapped  within  the  high  buoyancy  frequency  region  of 
the  inversion,  and  wave  propagation  occurred  in  nearly  all  horizontal 
directions.  Phase  speeds  ranged  from  less  than  1  m/sec  to  over  20  m/sec,  and, 
as  in  previous  studies,  no  evidence  for  critical  layers  was  found.  Stilke 
(1973)  observed  ducted  internal  gravity  waves  over  Hamburg,  Germany,  using  an 
instrumented  tower.  Wave  periods  ranged  from  approximately  4  to  20  minutes 
with  peak  to  peak  amplitudes  of  up  to  100  meters.  Sethuraman  (1977,  1980) 
observed  internal  waves  using  an  instrumented  tower  at  Long  Island,  New  York. 
A  strong  inversion  was  present  due  to  the  advection  of  warm  air  from  land  over 
the  cold  water.  These  waves  had  periods  of  order  5  to  10  minutes,  and 
Sethuraman  estimated  wavelengths  and  phase  speeds  to  be  in  the  4  km  and 
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4  m/sec  range.  The  phase  speeds  were  far  greater  than  the  mean  wind,  and  this 
excluded  the  existence  of  critical  layers.  Sethuraman  observed  periods  of 
persistent  high  levels  of  turbulence,  lasting  well  over  an  hour,  and  this  is 
in  contrast  to  the  frequently  observed  short-lived  turbulent  bursts.  The 
measurements  were  not  detailed  enough  to  yield  information  of  the  generation 
mechanism  or  the  physics  of  the  persistent  turbulence,  and  Sethuraman  did  not 
speculate  on  either  of  these. 


3.3  Characteristic  Features  of  Long  Internal  Gravity  Waves 

These  observations  demonstrate  that  gravity  waves  are  common  in  a  stably 
stratified  atmosphere.  Moreover,  many  of  the  waves  that  were  observed  were 
significantly  longer  than  height  of  the  tropopause  for  freely  propagating 
waves,  or  than  the  duct  height  for  ducted  waves.  In  this  section  we  discuss 
some  of  the  characteristic  features  of  long  internal  gravity  waves,  so  that 
they  may  be  identified  easily.  In  particular,  we  consider  the  following 
questions: 

(i)  Do  these  waves  have  any  distinguishing  features  (e.g.,  shape)  that 
would  help  to  identify  them? 

(ii)  How  long  do  they  persist? 

The  question  of  their  generation  has  already  been  discussed.  The  relation 
between  the  amplitude  of  an  internal  gravity  wave  and  the  intensity  of  the 
resulting  radar  signal  may  be  obtained  from  Appendix  A. 

The  equation  governing  the  horizontal  propagation  of  long  wavelength, 
small  amplitude,  two-dimensional  waves  ■>s  the  linear  wave  equation, 


ii±  .  a2 

3t2  3X2 


(3.12) 


where  a  is  the  speed  found  from  Eq.  (3.1)  for  the  appropriate  mode.  [Recall 
that  the  vertical  structure  was  determined  by  Eq.  (3.1).]  The  waves  are  a  few 
kilometers  long,  and  Eq.  (3.12)  governs  their  propagation  over  distances  of 


tens  of  wavelengths  (e.g.,  distances  less  than  100  km).  In  particular, 
Eq.  (3.12)  applies  over  the  field  of  vision  of  most  radars.  It  is  well  known 
that 


4>(x,t)  =  f(x-at)  (3.13) 

solves  Eq.  (3.12)  for  any  function  f.  Because  f  is  arbitrary  (being 
determined  by  initial  conditions),  it  follows  that  no  particular  wave  shape  is 
preferred  by  Eq.  (3.12).  Thus,  the  simplest  answer  to  question  (i),  above,  is 
negative. 

In  the  atmosphere,  long  internal  gravity  waves  that  are  approximately 
two-dimensional  have  been  observed  to  propagate  over  hundreds  of  kilometers 
(e.g.,  Christie,  Muirhead,  and  Hales,  1979;  Clarke,  Smith,  and  Reid,  1981). 
Equation  (3.12)  is  inadequate  to  describe  the  propagation  of  these  waves  over 
very  long  distances,  because  small  but  cumulative  effects  that  were  omitted 
from  Eq.  (3.12)  eventually  have  enough  time  to  exert  a  measurable  influence. 
For  example,  Grimshaw  (1980)  has  shown  that  in  a  moist,  compressible 
atmosphere,  long  waves  of  small  amplitude  propagate  according  to  Eq.  (3.12)  on 
a  short  time  scale,  and  according  to  the  Korteweg-deVries  equation. 


uT  +  6uux  +  uxxx  =  0  ,  (3.14) 

on  a  long  time-scale.  Here  x  =  x  -  at,  x  is  a  slow  time  (or  space)  scale,  and 
u  is  proportional  to  <j>  in  Eq.  (3.12).  The  assumptions  underlying  Eq.  (3.14) 
are  that  weak  nonlinearity  and  weak  dispersion  are  of  equal  importance,  and 
that  both  dominate  weak  three-dimensionality,  weak  dissipation  and  slow 
changes  in  topography  (including  the  height  of  the  tropopause). 

If  weak  three-dimensionality  is  as  important  as  the  nonlinearity,  then 
Eq.  (3.14)  should  be  replaced  by  the  Kadomtsev-Petviashvil i  equation. 
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(3.15) 


(u„  +  6uu  +  u  )  +  u  „  =  0  , 

'  t  uux  uxxx'x  nn  ’ 

where  n  is  measured  along  the  crest  of  the  wave  (e.g.,  see  Kadomtsev  and 
Petviashvil  i ,  1970;  or  Ablowitz  and  Segur,  1979).  Other  balances  are 

possible  as  well. 

Equations  (3.14)  and  (3.15)  are  among  a  class  of  special  nonlinear 
partial  differential  equations  that  can  be  solved  exactly.  In  particular, 
both  equations  have  N-soliton  solutions  [for  (3.14),  see  Gardner,  Greene, 
Kruskal ,  and  Miura,  1967,  1974;  for  (3.15),  see  Zakharov  and  Shabat,  1974; 
and  Satsuma,  1976],  The  one-sol iton  solution  for  Eq.  (3.14)  is  simply  a 
sol itary  wave: 


u(x»t)  =  sech2  j  K(x  -  K2  t  +  x0) 


(3.16) 


where  K  and  x0  are  arbitrary  parameters.  This  solution  is  shown  in 
Figure  3.6.  It  is  also  a  solution  of  Eq.  (3.15).  In  both  cases,  it 
represents  a  plane  wave,  traveling  without  change  of  form  in  the  x-direction. 
In  the  original  coordinate  system  (x,t),  it  travels  in  the  x-direction  with  a 
speed  slightly  faster  than  a  in  Eq.  (3.12),  and  travels  as  a  permanent  wave 
over  distances  much  longer  than  the  range  of  validity  of  Eq.  (3.12). 

The  relation  between  Eqs.  (3.12)  and  (3.14)  may  be  stated  as  follows. 
Consider  a  long,  two-dimensional  wave  of  small  amplitude.  Its  vertical 
structure  is  given  by  Eq.  (3.1),  but  its  horizontal  structure  may  be  specified 
arbitrarily,  provided  only  that  it  is  made  up  entirely  of  long  waves.  It  will 
propagate  as  a  wave  of  permanent  form,  according  to  Eq.  (3.12),  over 
relatively  short  distances  (less  than  100  km,  say).  Over  much  longer 
distances,  most  of  these  waves  will  evolve  slowly,  according  to  Eq.  (3.14)  if 
weak  nonlinearity  and  dispersion  dominate  in  other  small  effects.  The  only 
waves  that  retain  their  shape  over  these  large  distances  are  solitary  waves, 
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given  by  Eq.  (3.16).  A  spectacular  example  of  such  a  two-dimensional  solitary 
wave  in  the  atmosphere  is  discussed  by  Clarke,  Smith,  and  Reid  (1980).  [For  a 
detailed  experimental  study  of  the  accuracy  of  Eq.  (3.14)  in  the  closely 
related  problem  of  long  gravity  waves  in  shallow  water,  see  Hammack  and  Segur, 
1974,  1978.] 

Now  let  us  consider  weakly  three-dimensional  waves,  and  Eq.  (3.15).  A 
two-sol iton  solution  of  Eq.  (3.15)  is 


u(x»n,t) 


An  U  +  expUJ  +  exp(<p2)  +  exp(*1  +  <t>2  +  A12)} 


(3.17) 


where 


K,  (x  +  Pj  n  -  Ci  t)  ,  Ci  =  K|  +  p| 


i  =  1,  2 


and 


exp(Ai2) 


3(Ki  -K2)A  -  (Pl  -p2)2 
3(Ki  +  K2)2  -  (Pl  -p2)2 


The  pattern  of  the  wave  crests  in  the  x_n  plane  (i.e.,  in  the  x-y  plane, 
scaled  appropriately)  is  shown  in  Figure  3.7,  for  a  particular  choice  of  the 
parameters.  This  entire  pattern  moves  in  the  x-direction  as  a  wave  of 
permanent  form.  More  complicated  wave  patterns  have  been  considered  in  detail 
by  Anker  and  Freeman  (1978). 


We  are  not  aware  of  any  observations  in  the  atmosphere  of  long  internal 
gravity  waves  that  resemble  the  pattern  shown  in  Figure  3.7.  However, 
Eq.  (3.15)  also  describes  the  evolution  of  surface  gravity  waves  in  shallow 
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water  (here  "shallow  water"  is  equivalent  to  "long  waves").  In  shallow  water 
such  waves  have  been  observed,  as  may  be  seen  in  Figure  3.8.  They  occur 
commonly  on  beaches  of  very  gentle  slope,  with  waves  that  are  too  small  to 
break.  There  seems  to  be  no  reason  why  they  should  not  occur  in  the 

atmosphere  as  wel 1 . 

An  interesting  feature  of  the  waves  shown  in  Figure  3.8  is  that  these 
waves  are  breaking  only  in  the  relatively  small  region  of  intersection,  where 
the  wave  amplitude  is  larger  (roughly  because  the  amplitudes  of  the  two  waves 
add  in  this  region).  The  level  of  turbulence  in  this  localized  breaking 

region,  therefore,  should  be  significantly  higher  than  anywhere  else  in  the 
entire  wave  pattern.  Consequently,  a  measuring  device  that  responds  only  to 
very  short  wavelengths,  such  as  a  radar,  is  likely  to  ignore  the  long-crested 
waves,  and  to  "see"  only  the  localized  region  of  intersection.  Thus  the  radar 
return  would  show  a  three-dimensional  target  even  though  the  waves  that 
produced  it  are  primarily  two-dimensional. 

Generating  mechanisms  for  these  waves  were  discussed  in  the  previous 

section.  From  the  viewpoint  of  the  radar  operator,  however,  the  specific 

generating  mechanism  may  be  unimportant.  Because  the  waves  propagate  over 
long  distances,  they  need  not  be  related  to  any  features  within  the  operate; ' s 
field  of  vision,  such  as  anomolies  in  the  temperature  of  the  sea  surface.  He 
sees  only  waves  that  propagate  across  his  field  of  vision  without  any 
significant  change  of  form.  If  the  waves  were  generated  only  moderately  far 
away,  their  shape  will  bear  the  signature  of  the  generating  mechanism.  If 

they  were  generated  very  far  away,  they  may  have  evolved  into  solitons.  There 

seems  to  be  no  reason  to  expect  that  solitons  should  appear  more  often  than 

other  long  waves,  however,  unless  it  is  known  that  all  sources  of  these  waves 
are  very  distant.  [Conceivably  an  operator  in  Hawaii  would  see  solitons 

relatively  often.]  From  the  operator's  standpoint,  the  most  important 
features  of  these  long  internal  gravity  waves  are  that  they  travel  relative  to 
the  wind  at  all  levels,  and  that  they  evolve  very  slowly,  if  at  all. 
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(Photograph  courtesy  of  Toedtmeier  and  No  Nukes  Calendar,  1979. 


Turbulence  intensities  produced  by  such  waves  are  likely  to  be  associated 
with  the  incremental  shear  added  by  the  wave  forcing  the  local  value  of  Ri  to 
decrease  below  0.25  leading  to  the  K-H  instability  detailed  in  Appendix  A. 

3.4  Review:  Possible  Sources  of  Clear  Air  Radar  Returns 
in  the  Marine  Boundary  Layer 

The  preceeding  discussion  has  focused  on  the  various  mechanisms  which  can 
give  rise  to  clear  air  returns  and  details  some  observations  of  such 
phenomena.  We  conclude  this  section  with  a  review  of  these  phenomena  in 
outline  format.  The  phenomena  are  first  grouped  according  to  whether  they 
will  or  will  not  travel  with  the  ambient  wind.  This  is  followed  by  a  cross 
referencing  according  to  various  other  observable  features.  Our  goal  is  to 
provide  the  reader  with  a  way  to  discriminate  between  flow  types  based  on  the 
observed  characteristics  of  the  radar  target. 


I.  Clear  air  returns  which  travel  with  the  local  wind  velocity. 

A.  Convective  phenomena  associated  with  the  statically  unstable 

boundary  layer. 

(1)  Individual  convective  cells;  horizontal  size:  50  to  2000 
meters,  lifetime:  20  to  40  minutes. 

(2)  Benard  cell  pattern,  collective  motion  of  individual  cells, 
size:  5  to  10  km;  lifetime:  several  hours. 

(3)  Thermal  streets,  collective  motion  of  individual  convective 

cells,  bands  of  thermals  aligned  with  mean  wind;  spacing 

between  bands:  1  to  8  km;  lifetime:  several  hours. 

B.  Kel vin-Helmholtz  Instabilities. 

(1)  Unstable  waves  associated  with  mean  shears  present  in  the 
atmospheric  boundary  layer.  Wavelengths,  periods  and  phase 
velocities  of  the  order  100  to  400  meters,  100  to  300  seconds, 
and  1  to  5  m/sec,  respectively,  have  been  observed. 


II.  Clear  air  returns  which  travel  at  velocities  different  from  the  local 
wind.  (Wave  phenomena  which  can  result  in  localized,  propagating 
regions  of  turbulence.) 

A.  Tropospheric  waves,  propagating  internal  gravity  waves  which  extend 
throughout  the  troposphere.  Observed  horizontal  wavelengths;  5  to 
400  km;  periods  5  to  100  minutes;  and  phase  speeds:  5  to  50 
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m/sec.  Common  mechanisms  for  wave  excitation  are  storms  and 
fronts. 


B.  Ducted  waves.  Internal  waves  which  are  confined  to  the  boundary 
layer  due  to  strong  inversions.  Wavelengths:  50  to  1000  meters, 
periods:  30  to  500  seconds,  and  phase  speeds:  1  to  20  m/sec. 

C.  Solitons.  Long  wavelength  internal  gravity  waves  of  moderately 
large  amplitude  which  are  shape-preserving  over  long  distances. 

D.  Individual  convective  cells  in  the  presence  of  strong  vertical  wind 
shear. 


III.  Radar  Returns  which  are  Localized 

A.  Individual  convective  cells.  Scale  of  such  thermals  ranges  from  50 
to  2000  meters  depending  on  position  within  boundary  layer. 

B.  Localized  turbulence  associated  with  internal  gravity  waves. 

Regions  of  small  scale  shear  instabilities  can  appear  which  are 
associated  with  propagating  waves.  Computations  show  that  the  size 
of  turbulence  region  should  be  of  order  10  (A/Au)  CBhase>  where  1 
and  Au  are  the  thickness  and  velocity  difference  or  the  locally 
unstable  region,  and  Cp^  is  the  phase  velocity.  Typical  values  for 
these  parameters  yield  sizes  of  order  100  to  500  meters. 

C.  Solitons.  The  possible  intersection  of  two  solitary  waves  could 
lead  to  localized  breaking  and  large  Cp. 

IV.  Radar  returns  which  are  relatively  short-lived. 

A.  Individual  convective  cells.  Lifetimes  range  from  20  to  40 
minutes,  but  near  the  end  of  its  existence,  cell  begins  to  break 
up.  Therefore,  one  may  only  have  a  trackable  entity  for  times  of 
order  20  minutes. 

B.  Localized  turbulence  associated  with  internal  waves,  same  as 
category  III.B.  Turbulence  lifetimes  are  of  order  10  aa/au,  and 
typical  values  yield  times  of  order  1  to  10  minutes. 


V.  Radar  returns  which  can  be  observed  through  the  well-mixed  boundary 
layer. 

A.  Individual  convective  cells  are  observable  due  to  increased  at 
their  boundaries,  and  may  be  tracked  throughout  the  boundary  layer. 
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VI.  Radar  returns  which  will  be  most  prominent  near  inversions,  or 
mixed- layer  tops. 

A.  All  types  due  to  the  increased  humidity  gradients  typically  found 
in  this  region. 


VII.  Radar  returns  influenced  by  surface  conditions. 

A.  Convective  activity  depends  strongly  on  surface  heating.  Any 
surface  mechanism  which  creates  localized  "hot  spots"  should  lead 
to  individual  convective  cells,  while  the  creation  of  extended 
lengths  of  warm  water  over  the  ocean  could  give  rise  to  collective 
motion  of  thermals. 

B.  Topographic  features  or  convective  activity  may  trigger  waves. 
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IV.  Preliminary  Analysis  of  NRL  Results  on  in  Situ  Measurements 
of  Some  Clear  Air  Returns 


Radar  returns  from  atmospheric  disturbances  were  recorded  during  a  series 
of  field  tests  conducted  at  the  AUTEC  Weapons  Range  during  May  13  to  May  17, 
1979.  According  to  the  Data  Catalogue  of  Aircraft  Measurements,  prepared  by 
Dr.  Lothar  Ruhnke  of  NRL,  radar  returns  were  observed  during  Flight  #7 
(afternoon  of  May  16)  and  Flight  #8  (morning  of  May  17).  As  seen  from 
Chapter  III  the  characteristic  most  likely  to  provide  discrimination  between 
targets  generated  by  wave  breaking  and  those  generated  by  convection  is  the 
speed  of  the  target  relative  to  the  wind.  In  this  chapter,  we  estimate  the 
difference  between  convectively  generated  targets  and  wave  generated  targets 
that  may  have  been  present  during  these  two  flights.  We  will  show  that  there 
exist  neutrally  stable  gravity  waves  in  the  atmosphere,  the  characteristics  of 
which  are  consistent  with  the  observed  radar  returns,  and  which  travel  at  the 
reported  speeds.  We  emphasize,  however,  that  the  data  available  are  sparse, 
so  that  these  calculations  are  only  suggestive. 

Figures  (4.1a,b)  and  (4.2a,b)  show  the  vertical  distributions  of  virtual 
potential  temperature  and  specific  humidity  that  were  measured  in  the  lower 
atmosphere  during  Flights  #7  and  #8,  respectively.  [The  logarithmic 
derivative  of  the  virtual  potential  temperature  is  approximately  t$(z)  in 
Eq.  (3.9).]  Note  that  in  both  cases  there  is  the  suggestion  of  a  thin 
unstable  layer  next  to  the  surface  (as  much  as  a  few  10's  of  meters  thick) 
capped  by  a  mixed  layer  (i.e.,  ti  «  0)  in  the  lower  800  meters  or  so  of  the 
atmosphere,  with  a  stably  stratified  (^  >  0)  region  above  it.  During 
Flight  #7  the  "mixed"  layer  is  slightly  stable  but  considerably  less  stable 
than  the  region  above  it.  The  stability  measured  above  the  mixed  layer 
(~  3°C/km)  is  comparable  to  typical  values  in  the  entire  troposphere.  (For 
example,  the  lapse  rate  for  the  standard  atmosphere  is  about  3.5°C/km.)  The 
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Figure  4.1(a)  -  Vertical  distribution  of  virtual  potential  temperature  durinq 
Flight  #7.  Taken  from  Ruhnke  (1979). 
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inversions  capping  the  mixed  layer  on  these  two  days  should  be  considered 
weak,  in  comparison  with  much  stronger  inversions  that  have  been  measured 
(e.g.,  Metcalf,  1975). 


As  discussed  in  Chapter  III,  vf  there  were  a  very  strong  inversion  at  the 
top  of  the  mixed  layer,  then 

(i)  the  gravity  waves  would  be  ducted  within  this  layer,  with 
exponentially  decreasing  amplitudes  above  it;  and 

(ii)  their  speeds  would  be  determined  primarily  by  the  distribution  of 
virtual  potential  temperature  in  the  lowest  1  or  2  kilometers  of  the 
atmosphere  (i.e.,  by  the  measured  distributions). 

However,  because  the  measured  inversions  are  so  weak,  neither  of  these 
statements  is  true.  In  particular,  the  measured  distributions  of  stability  in 
Figures  (4.1a)  and  (4.2a)  are  inadequate  to  estimate  wave  speeds.  The  wave 
speeds  we  calculate  below  depend  on  an  assumed  distribution  of  stability  in 
the  entire  troposphere  and  on  an  assumed  location  of  the  tropopause. 

The  observed  features  of  the  radar  targets  were  the  following: 

(i)  Location:  They  were  observed  near  the  top  of  the  mixed  layer  (about 
800  m) . 

(ii)  Size  and  shape:  The  horizontal  extent  of  a  disturbance  was  5  km  or 
more.  (An  aircraft  flying  at  about  50  m/sec  traversed  the  target  in 
about  100  seconds.)  Some  of  the  returns  had  a  "crescent"  shape, 

(iii)  Speed:  Some  of  the  targets  travelled  at  speeds  different  from  the 
winds.  The  available  information  on  wind  speed  and  direction  in  the 
lowest  2  km  of  the  atmosphere  is  shown  in  Table  4.1.  Typical  wind 
speeds  did  not  exceed  about  6  m/sec  (~  12  knots).  Some  radar 
targets,  on  the  other  hand,  travelled  at  speeds  of  10-15  m/sec  and 
in  directions  significantly  different  from  the  wind.  Thus  these 
targets  travelled  at  speeds  on  the  order  of  5-10  m/sec 
relative  to  the  winds. 
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Table  4.1  -  Vertical  distribution  of  wind  direction  (degrees  from  North) 
and  wind  speeds  (knots)  measured  at  AUTEC  during  May  14-17,  1979  Data 
from  D.  F.  Hemenway. 
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There  is  a  high  probability  that  numerous  convective  cells  existed  on 
these  two  occasions.  Figures  4.1a  and  4.2a  show  that  an  anomalous  warm  spot 
on  the  surface  of  as  much  as  a  few  tenths  of  a  degree  should  be  adequate  to 
induce  a  buoyant  cell  with  sufficient  buoyancy  to  rise  to  heights  of  800 
meters  or  more.  The  evaporative  cooling  of  the  ocean  surface  suggested  by 
Figures  4.1b  and  4.2b  strongly  suggests  that  the  upper-most  layer  of  the  ocean 
is  also  unstably  stratified.  As  this  layer  randomly  turns  over  it  should 
occassional ly  provide  sufficiently  warm  surface  spots  to  induce  distinct 
convective  cells  in  the  atmosphere. 

A  well  organized  convective  plume  that  carries  surface  air  with  a 
humidity  of  approximately  0.016  gm/gm  to  an  elevation  of  800  meters  where  the 
ambient  humidity  is  approximately  0.013  gm/gm  may  be  expected  to  produce  C(,  in 
the  outer  shell  of  the  cell  at  this  altitude  of  as  much  as  10-12,  as  one  may 
infer  from  Eq.  (2.7).  Dilution  of  the  surface  moisture  as  the  cell  rises 
through  the  mixed  layer  may  be  expected  to  reduce  this  value,  but  there 
appears  to  be  ample  potential  for  it  to  reach  detectable  levels. 

The  wind  data  provided  in  Table  4.1  indicates  that  there  was  little 
vertical  wind  shear  below  1  km,  so  that  any  convective  cells  reaching  the 
level  of  800  meters  should  be  expected  to  move  with  the  ambient  wind  speed 
regardless  of  how  the  convective  cell  is  generated  at  the  surface.  Thus  the 
observations  are  completely  consistent  with  convectively  generated  targets, 
except  for  the  report  that  the  targets  moved  at  speeds  and  direction 
distinctly  different  from  that  of  the  wind.  This  suggests  the  possibility  of 
wave-generated  targets. 

Because  of  the  lack  of  any  atmospheric  data  above  2  km,  we  will  use  a 
very  simplified  model  of  the  atmosphere  to  calculate  wave  speeds.  Our 
assumptions  are  these. 

(i)  The  height  of  the  tropopause  is  10  km,  and  acts  like  a  rigid  lid; 
i.e.,  F  *  0  at  z  *  10  km  for  Eq.  (3.1). 

(ii)  The  sea  surface  acts  like  a  rigid  floor;  F  =  0  at  z  =  0. 

(iii)  The  stability  parameter,  B(z)  in  Eqs.  (3.1)  and  (3.9)  is  given  by: 


54 


0 


0  <2  <  800  m 


(4.1) 


‘  = 

I  1.2  x  10-4  m-1  ,  800m  <  z  <  10  km 

(iv)  The  wind  speed  and  direction  are  constant.  Consequently,  we  may 
rewrite  Eq.  (3.1)  in  a  coordinate  system  traveling  with  the  wind,  so 

that  the  computed  wave  speed  will  be  the  speed  relative  to  the  wind, 

(v)  We  may  use  the  Boussinesq  approximation  [p0  =  constant  in 
Eq.  (3.2)]. 

Probably  none  of  these  assumptions  precisely  simulates  the  atmospheric 

conditions  on  May  16  or  May  17,  1979.  However,  they  should  be  sufficiently 

accurate  that  the  resulting  wave  speeds  may  be  considered  representative. 

With  these  assumptions,  Eq.  (3.1)  becomes: 


with  e(z)  given  by  Eq.  (4.1).  Because  the  observed  wavelengths  exceed  5  km, 
we  will  also  take  the  long  wave  limit  (k  0).  Now  Eq.  (4.2)  may  be  solved  by 
elementary  means.  It  has  an  infinite  nunber  of  solutions,  each  representing  a 
different  vertical  mode.  Table  4.2  gives  the  vertical  structure  and  wave 
speed  (relative  to  the  wind)  of  the  first  seven  vertical  modes.  We  note 
several  points. 

(i)  If  the  rigid  lid  at  the  tropopause  were  replaced  with  a  free 
surface,  the  results  in  Table  4.2  would  not  change  significantly. 
An  additional  mode  (#0)  would  enter  at  the  top  of  the  table,  with  a 
speed  significantly  larger  than  35  m/sec. 

(ii)  The  first  three  or  four  modes  would  be  virtually  unaffected  if  the 
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Table  4. 


Vertical 

mode 

number 


P 


-  Vertical  structure  and  speeds  of  the  lowest  modes,  for  Eq.  (4.2). 


Distribution  of 
vertical  velocity 
below  the  tropopause 
(at  a  fixed  time  and 
horizontal  position) 


Height  of 
lowest  ■*  Wave 

maximum  in  speed 

horizontal  relative 

Distribution  of  shear,  to  the 

Horizontal  velocity  l3u/9z|  wind 

in  the  lowest  2  km  (m)  (m/sec) 


lO.OOOm 


z  *  2000  m 
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mixed  layer  below  800  m  were  removed. 

(iii)  Modes  4-7  all  have  wave  speeds  in  the  range  of  those  observed 
(roughly  5-10  m/sec  relative  to  the  wind). 

(iv)  For  all  modes  beyond  about  5,  the  long  wave  approximation  holds  for 
any  wave  significantly  longer  than  about  1  or  2  km.  Waves  in  these 
modes  with  horizontal  wavelengths  on  the  order  of  5  km  or  more 
certainly  would  propagate  like  long  waves. 

(v)  Modes  5-7  (and  all  higher  modes)  possess  regions  of  maximal  shear 
within  a  few  hundred  meters  of  the  top  of  the  mixed  layer.  Such 
regions  are  capable  of  generating  radar  returns,  by  the  scenario 
given  in  the  following  paragraph. 

The  data  in  Table  4.1  are  consistent  with  vertical  wind  shears  of  the 
order  of  1  to  2  knots  over  a  100  ft  (30  meters)  altitude  change  at  the  top  of 
the  mixed  layer,  where  Figures  4.1a  or  4.2a  show  the  potential  temperature 
change  to  be  of  the  order  of  0.1  to  0.3°C.  These  estimates  yield  a  Richardson 
number  of  the  order  of  0.1  to  1.0.  It  is  thus  quite  possible  for  there  to  be 
a  region  where  Ri  «  0.3.  The  addition  of  the  shear  associated  with  wave  modes 
5  to  7  with  even  modest  amplitudes  would  then  be  sufficient  to  reduce  the  Ri 
below  0.25  and  induce  the  type  of  wave  breaking  calculated  in  Appendix  A. 
This  region  of  breaking  wave  would  travel  with  the  wave  as  long  as  the  time 
scale  of  the  wave  breaking  is  less  than  the  time  required  for  the  wave  to 
pass.  Figure  A. 10  shows  that  the  time  scale  over  which  there  appears  a  large 
increase  in  temperature  variance  (which  should  be  the  same  as  that  for 
humidity  variance)  is  approximately  10  A/Au.  For  a  shear  layer  30  meters 
(100  feet)  thick  with  a  Au  of  1  m/sec  (2  knots)  this  gives  a  lifetime  of 
approximately  300  sec.  This  would  provide  a  region  1.5  to  3.0  km  long  of 
enhanced  Cp  travelling  with  a  long  wave  of  speed  5  to  10  m/sec  relative  to  the 
wind.  From  Figure  A. 9  the  peak  value  of  H ’ H '  should  be  approximately  0.3  AH2. 
The  moisture  increments  in  Figures  4.1b  and  4.2b  just  above  800  m  would  then 
be  sufficient  to  yield  values  of  H ' H *  «  10-6  which  in  turn  through  Eq.  (2.7) 
would  lead  to  ~  10-13,  a  borderline  value  for  detectability. 


57 


These  calculations  demonstrate  the  possibility  that  long  internal  gravity 
waves  were  present  on  May  16  and  17,  1979,  that  they  travelled  relative  to  the 
wind  at  speeds  of  5  to  10  m/sec,  that  they  generated  regions  of  high  (i.e., 
radar  targets)  near  the  top  of  the  mixed  layer,  and  that  these  targets  also 
travelled  relative  to  the  wind.  As  discussed  earlier,  targets  generated  by 
buoyant  plumes  also  were  possible  on  these  days.  Because  the  inversion 
measured  on  these  two  days  was  fairly  weak,  one  would  expect  typical 
plume-generated  targets  (travelling  with  the  wind)  to  be  somewhat  stronger 
than  typical  wave-generated  targets  (travelling  different  from  the  wind).  If 
the  inversion  were  stronger,  this  discrepancy  in  strength  of  signal  would  be 
mitigated. 


i 
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V.  Concluding  Remarks 


We  have  provided  a  catalogue  of  the  dominant  types  of  fluid  dynamic 
features  which  may  be  expected  to  yield  trackable  radar  targets  in  the 
atmospheric  boundary  layer.  The  simplest  distinguishing  feature  between 
convectively  produced  targets  and  targets  produced  by  wave  breaking  is  the 
velocity  at  which  the  target  travels  relative  to  the  wind.  Except  in  the 

presence  of  strong  vertical  wind  shear,  convectively  generated  targets  will 
travel  with  the  ambient  wind,  while  wave  generated  targets  may  or  may  not 
travel  with  the  wind.  Conditions  favorable  to  convectively  generated  targets 
occur  quite  frequently;  those  for  wave  generated  targets  are  likely  to  occur 
somewhat  less  frequently. 

We  have  analyzed  in  detail  one  particularly  interesting  type  of  wave 
generated  target  which  travels  at  velocities  significantly  different  than  the 
ambient  wind.  A  crucial  part  of  this  scenario  is  the  Kel vin-Helmholtz  wave 

breaking  which  occurs  in  a  fluid  shear  layer  where  the  Richardson  nunber  drops 
below  0.25.  This  process  is  studied  in  detail  in  Appendix  A  to  show  the  life 
history  of  temperature  and  humidity  fluctuations  in  such  a  wave-breaking 
process.  The  turbulent  transport  model  used  allows  the  buildup  and  subsequent 
decay  of  these  variances  to  be  simulated  as  a  function  of  the  initial 
conditions  of  the  shear  layer.  The  results  of  this  detailed  numerical 

simulation  were  used  in  our  analysis  of  what  may  have  been  occurring  during 

the  NRL  tests. 

This  initial  study  has  been  primarily  involved  with  the  general  features 
of  clear  air  radar  targets.  We  recommend  a  second  phase  to  investigate 
particular  ways  in  which  these  signals  may  become  amplified.  Interesting 
questions  which  should  be  answerable  by  further  analysis  are: 

1)  Are  there  conditions  for  which  individual  radar  targets  can  be 
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traced  back  to  the  atmospheric  disturbance  which  induced  it  in  a 
coherent  fashion? 

2)  The  complimentary  question:  Under  what  conditions  can  a  discrete 
atmospheric  disturbance  be  amplified  sufficiently  to  produce  a 
detectable  target? 


APPENDIX  A 


Detailed  Examination  of  Kel vin-Helmhol tz  Wave  Breaking* 


A.l  Introduction 

Shear  instabilities  are  now  widely  recognized  as  a  significant  mechanism 
for  producing  turbulence  and  mixing  in  stably-stratified  fluids.  There  are 
numerous  observations  of  "billow  turbulence"  in  the  atmosphere  and  oceans  (see 
Maxworthy  and  Browand,  1975)  and  these  events  have  been  associated  with 
Kel vin-Helmholtz  instabilities  of  the  wind  profile.  Laboratory  studies  by 
Thorpe  (1973)  show  similar  vortex  roll-up  features  and  the  generation  of  a 
turbulent  layer  of  fluid  which  mixes  and  spreads  the  initial  shear  layer 
before  decaying  back  to  a  quiescent  flow. 

The  initial  stages  of  the  instability  are  now  well  understood.  Finite 
amplitude  numerical  computations  by  Patnaik,  et  al .  (1976),  Peltier,  et 
al .  (1978)  have  confirmed  the  linear  stability  predictions  of  growth  rates  and 
mode  structure  of  the  growing  waves,  and  have  gone  on  to  calculate  the  roll-up 
of  the  vortex  layer.  However,  these  laminar  calculations  at  moderate  Reynolds 
numbers  have  been  unable  to  identify  the  secondary  instability  which  results 
in  the  breakdown  of  the  vortex  into  a  turbulent  layer.  Davis  and  Peltier 
(1979)  made  calculations  up  to  a  Reynolds  number  of  500  based  on  shear  layer 
thickness  but  failed  to  produce  a  secondary  instability.  In  the  latter  paper, 
the  authors  calculated  local  Rayleigh  numbers  in  the  flow,  and  showed  that 
with  a  Reynolds  number  of  500,  a  region  of  very  high  Rayleigh  number  was 
produced  in  the  rolled-up  vortex.  They  speculate  that  this  region  is  in  fact 
convectively  unstable,  but  that  the  most  unstable  models  will  be  longitudinal 


*  This  work  was  also  partially  supported  by  the  Naval  Air  Systems  Command 
under  Contract  No.  N00019-80-C-0617. 
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rolls  since  the  flow  is  strongly  sheared.  This  hypothesis  accounts  for  the 
fact  that  the  instabilities  have  not  been  triggered  in  the  two-dimensional 
numerical  models. 

The  aim  of  the  present  study  is  to  model  the  breakdown  of  the  rolled-up 
vortex  layer  and  the  subsequent  spread  and  decay  of  the  turbulence  using  a 
second-order  closure  scheme  to  describe  the  three-dimensional,  small  scale 
turbulence  field.  Second-order  closure  models  have  been  used  extensively  to 
study  hori zontal ly-homogeneous  atmospheric  boundary  layer  flows,  and  in 
particular  are  capable  of  representing  the  major  dynamical  processes  involved 
in  turbulent  convection  (see  Wyngaard  and  Cote,  1974;  Lewellen,  Teske  and 
Donaldson  1976;  Yamada  and  Mellor,  1975).  [There  are  difficulties  associated 
with  the  turbulent  transport  of  turbulent  kinetic  energy  (see  Zeman  and 
Lumley,  1976;  Andre,  et  al.,  1976),  but  this  is  not  a  serious  problem  for  the 
prediction  of  the  mean  profiles  of  the  major  second-order  quantities.]  Since 
the  model  calculates  the  growth  of  turbulent  correlations  in  an  unstable 
environment  moderately  accurately  (because  the  growth  is  driven  directly  by 
production  terms  in  the  equations  without  recourse  to  empirical  closure),  we 
may  reasonably  expect  the  results  of  the  calculation  to  have  some  validity. 
One  would  not  have  the  same  expectation  of,  say,  a  mixing-length  closure 
model,  since  the  time  scale  of  the  growth  of  the  turbulence  is  the  same  as  the 
mean  flow  time  scale,  while  the  mixing-length  model  implicitly  assumes  that 
the  turbulence  is  in  local  equilibrium.  Thus  turbulent  stresses  and  heat 
fluxes  would  be  generated  as  soon  as  the  density  profile  became  unstable  if 
the  mixing-length  were  related  to  the  stability  as  in  the  model  of  Orlanski 
and  Ross  (1973). 

As  we  shall  demonstrate,  the  second-order  closure  model  produces  results 
which  are  certainly  qualitatively  correct.  Quantitative  data  on  turbulence 
quantities  is  not  sufficiently  detailed  to  make  any  meaningful  assessment  of 
model  accuracy,  but  there  are  more  general  quantities  associated  with  the 
breakdown  which  can  be  used  for  comparison.  There  are  estimates  of  the  time 
scale  of  the  turbulence  decay  and  measurements  of  the  final  states  from  the 
experiments  by  Thorpe  (1973),  for  example.  The  final  state  of  the  shear  layer 
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is  an  important  quantity  in  deciding  the  importance  of  Kel vin-Helmholtz 
instability  as  a  mixing  process  in  stratified  fluids,  since  this  is  a  measure 
of  the  total  amount  of  mixing  produced  by  the  event.  Another  source  of  data 
is  high-resolution  radar  probing  of  the  atmosphere,  where  vortex  roll-up 
events  have  b«»en  reported  by  many  investigators,  as  discussed  in  Chapter  III 
of  the  main  report.  However,  although  it  is  possible  to  estimate  the 
intensity  of  the  turbulent  fluctuations  from  the  intensity  of  the  radar  echo, 
the  measurements  of  the  initial  wind  and  temperature  profiles  immediately 
prior  to  the  instability  are  generally  inadequate  for  a  crucial  test. 

In  spite  of  the  lack  of  detailed  information,  there  is  certainly 
sufficient  data  to  determine  whether  the  turbulence  model  is  reasonably 
accurate,  and  hopefully  some  of  the  numerical  results  will  prompt  more 
detailed  measurements. 


A. 2  Model  Equations  and  Numerical  Solution 

The  second-order  closure  scheme  employed  in  this  study  has  been  described 
in  detail  by  Lewellen  (1977)  and  we  therefore  present  only  the  final 
equations.  We  consider  incompressible,  Boussinesq  flow,  and  write  the 
Reynolds-averaged  equations  with  turbulent  correlations  denoted  by  an  overbar. 
For  two-dimensional  motion  in  the  (x,z)-plane,  with  z  in  the  upward  vertical 
direction  we  define  a  stream-function  v,  such  that  u  =  ay/az,  W  =  -  av/ax. 
Then  the  vorticity  4  =  v^y,  and  the  equations  of  motion  are 
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In  the  above,  (U,W)  is  the  mean  velocity,  T  is  the  mean  temperature, 
(u,v,w)  is  the  turbulent  velocity  fluctuation,  o  is  the  temperature 
fluctuation,  and  q  =  [u7  +  "v7  +  w7] A  is  a  length  scale,  and  appears  in 
the  modelled  terms  in  the  Reynolds-averaged  turbulent  correlation  equation. 
The  values  of  the  empirical  constants  are  as  follows:  b  =  1/d,  A  =  0.7b, 
Vc  =  0.3,  and  s  =  1.8.  These  values  were  previously  chosen  to  match  model 
results  with  a  set  of  experiments  with  simple  geometries  as  described  by 
lewellen  (1977). 

For  the  Kel vin-Helmholtz  instability  calculation,  we  begin  with  a  profile 
of  temperature  and  velocity,  viz, 
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tanh  — 
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(A.  11 ) 


U  =  all  tanh  - 
6 


(A. 12) 
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Note  that  velocity  and  temperature  both  have  the  same  profile  shape.  The 
initial  conditions  define  a  Richardson  number  Ri  =  ( g/T0 ) (&T6/aUz) ,  which  will 
be  the  basic  dimensionless  number  in  our  experiments.  The  problem  also 
contains  a  Reynolds  number,  but  since  we  have  omitted  explicit  laminar 
diffusive  terms  from  the  equations  of  motion,  we  are  effectively  studying  the 
high  Reynolds  number  limit. 

In  order  to  begin  the  integration,  we  have  to  specify  initial  values  for 
the  turbulence  energy,  qz,  and  the  length  scale  A;  in  the  calculations 
reported  here,  all  other  turbulence  quantities  were  initially  set  to  zero.  A 
perturbation  vorticity  was  also  added  so  that  the  instaDil ity  would  amplify 
this  disturbance  and  produce  the  vortex  roll-up.  The  specification  of  these 
initial  conditions  will  be  discussed  in  the  next  section. 

Periodic  boundary  conditions  were  employed  in  the  x-direction,  ana  the 
length  of  the  domain  was  chosen  somewhat  larger  than  the  wavelength  of  the 
fastest  growing  linear  mode.  This  choice  is  suggested  by  the  nonlinear 
calculations  of  Patnaik,  et  al.  (1976),  which  show  that  the  longer  waves  grow 
to  larger  amplitude,  and  may  therefore  dominate  the  finite  amplitude 
development  of  the  wave.  For  most  of  our  integrations,  a  domain  length  of  156 
was  used,  although  some  runs  were  made  with  different  lengths  to  investigate 
the  sensitivity. 

The  boundaries  in  the  z-direction  were  placed  sufficiently  far  from  the 
shear  layer  for  them  to  have  negligible  effect  on  the  flow  development.  A 
zero  normal  gradient  condition  was  specified  on  all  flow  quantities  on  the 
upper  and  lower  boundaries. 

The  equations  of  motion  are  discretized  on  a  finite-difference  grid  which 
is  uniformly  spaced  in  the  horizontal  but  has  finer  vertical  resolution  in  the 
vicinity  of  the  shear  layer.  The  equations  are  solved  using  centered  spatial 
derivatives,  and  an  ADI  scheme  for  the  temporal  derivative.  The  advection 
operator  conserves  first  and  second  moments  apart  from  time-dependence  errors 
arising  from  the  use  of  explicit  advective  velocity  components  in  the  ADI 
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scheme.  The  time-stepping  is  made  more  stable  by  the  introduction  of  a 
coupling  between  the  mean  variables  and  the  turbulent  fluxes  by  the  artifice 
of  an  additional  implicit  eddy  diffusive  flux  which  is  explicitly  subtracted 
when  the  correct  turbulent  flux  is  added  (lewellen  and  Sheng,  1980).  Provided 
the  time-step  is  sufficiently  small,  the  errors  introduced  by  this  procedure 
will  be  negligible. 

Finally,  the  Poisson  equation  for  the  streamfunction  is  solved  directly 
using  the  decomposition  method  of  Swarztrauber  and  Sweet  (1975). 

A. 3  Results 

We  consider  an  initial  state,  as  described  in  Section  A. 2,  given  by 


T  =  T0  +  aT  tanh  -  (A. 13a) 

U  =  aU  tanh  —  (A. 13b) 
6 

which  defines  a  Richardson  number  Ri  =  (gAAT/TQAuO  •  This  is  in  fact  the 
minimum  gradient  Richardson  nunber,  defined  by  Rig  =  £(g/T0) (3T/az)J/(dU/azj2, 
which  occurs  at  z  =  0.  Linear  inviscid  theory  (Miles,  1961)  shows  that  the 
flow  is  stable  to  small  perturbations  provided  Ri  >  0.85. 

We  first  present  detailed  results  from  a  calculation  with  Ri  =  0.1.  In 
this  example,  we  applied  a  sinusoidal  (in  x)  variation  to  the  vorticity  in 
order  to  initiate  the  instability.  The  perturbation  vorticity  amplitude  of 
0.1^u/6  is  quite  large;  this  is  because  we  are  mainly  interested  in  the 
finite-amplitude  breaking  of  the  billow  and  the  subsequent  turbulent 
development,  so  that  we  do  not  wish  to  spend  time  computing  the  initial  growth 
from  a  very  small  perturbation.  Checks  were  made  to  ensure  that  the  initial 
perturbation  was  sufficiently  small  that  it  did  not  influence  the 
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finite-amplitude  development. 


The  turbulence  energy  q2  was  initially  set  at  10-iJ(£U)^  throughout  the 
domain,  with  the  turbulence  length  scale,  A,  set  equal  to  0.46.  The 
turbulence  energy  is  initially  small  enough  to  be  effectively  zero,  but  the 
initial  length  scale  does  have  some  effect  on  the  results.  Variations  in 
these  initial  values  will  be  discussed  in  more  detail  in  Section  6  below. 

The  integration  was  performed  using  a  domain  of  length  15o  and  a  height 

of  126  centered  on  the  mid-point  of  the  shear  layer.  The  computational  mesh 

consisted  of  41  x  61  points,  with  41  points  spaced  uniformly  in  the 
horizontal,  and  a  non-uniform  vertical  mesh  giving  a  grid  spacing  of  O.lo  in 
the  central  region  and  expanding  out  to  roughly  U.2o. 

Figure  A.l  shows  a  sequence  of  contour  plots  of  the  dimensionless 
temperature  field,  (T  -  T0)/aT,  at  dimensionless  time  x  =  1.5,  2.8,  4.3,  5.8, 
7.3,  and  11.8,  where  x  =  t(gaT/T0£u).  Corresponding  plots  of  the  small-scale 
turbulent  energy  q ^/ ( aU ) ^  are  presented  in  Figure  A. 2.  The  bil low  develops 
initially  as  an  effectively  laminar  flow  with  a  growing  vortex  core  in  the 
center,  and  thin  braids  forming  between  the  cores.  The  initial  turbulent 

perturbation  in  q2  is  amplified  slowly  in  the  braids  where  the  shear  is 
highest,  and  the  local  Richardson  number  is  lowest.  At  x  =  2.8  the 

instability  has  vertically  reached  its  maximum  amplitude,  with  temperature 
contours  quite  convoluted  in  the  vortex  core,  and  the  braids  have  become  very 
thin.  At  this  stage  the  temperature  field  is  statically  unstable  in  the  core, 
and  the  q2  contours  show  growth  of  turbulence  associated  with  the  buoyant 
instability,  although  the  turbulence  is  still  relatively  small  everywhere  with 
the  buoyantly  produced  turbulence  in  the  core  comparable  with  the  turbulence 
levels  in  the  braids. 

At  the  next  time,  x  =  4.3,  most  of  the  temperature  structure  in  the  core 
has  been  mixed  by  the  turbulence  which  has  reached  a  relatively  high  level 
around  the  outer  part  of  the  core;  the  maximum  value  of  q2  is  0.18(MJ)2.  The 
core  has  already  begun  to  spread  horizontally,  although  the  braids  are  still 


69 


leths  of  dimensionless  temperature  (T  -  T0)/aT  for  the  case 
he  contour  interval  is  0.2,  and  positive  contours  (i.e., 
e  denoted  by  a  dashed  line,  (a)  r  =  1.5. 


igure  A.l (c) 


gure  A. 1(d) 
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Figure  A.l (e) 


Figure  A.1 (f ) 


gure  A. 2  -  Isopleths  of  dimensionless  turbulence  kinetic  energy,  q2/AU2  at 
e  same  times  as  Figure  A.l,  i.e.,  (a)  t  =  1.5,  contour  interval  is  0.001. 
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interval  is  0.01 


igure  A. 2(e)  -  t  =  7.3,  contour  interval  is  0.0007. 


interval  is  0.002. 


quite  distinct.  At  t  =  5.8,  the  braids  have  been  mixed,  and  the  turbulence 
has  spread  throughout  the  layer  in  the  horizontal.  There  is  evidence  of 
smaller  scale  disturbances  in  the  temperature  field,  and  a  movie  of  the  time 
evolution  shows  them  to  be  travelling  wave-like  disturbances  moving  along  the 
top  and  bottom  of  the  layer,  and  also  weakly  rolling  up  in  a  manner  similar  to 
the  original  instability.  However,  they  are  also  being  mixed  by  the  high 
background  turbulence  so  that  they  do  not  persist  long.  There  is  also  some 
remnant  of  the  original  large  scale  vortex  which  continues  to  slowly  turn  the 
isotherms.  We  can  still  identify  the  regions  of  production  of  the  turbulent 
energy,  since  the  level  of  q2  in  the  center  is  about  twice  the  level  at  the 
edges  of  the*  domain. 

The  trend  toward  flattening  the  isotherms,  and  spreading  the  turbulence 
horizontally  continues  until  at  the  latest  time,  t  =  11.8,  there  is  virtually 
no  further  mean  motion,  and  the  turbulence  is  almost  homogeneous  in  the 
horizontal.  The  turbulence  level  has  dropped  to  0.022(aU)2  at  this  stage, 
i.e.,  almost  a  factor  of  10  below  its  maximum  amplitude  which  occurred  around 
i  =  4.3. 

The  evolution  described  above  is  qualitatively  similar  to  the  turbulent 
breakdown  of  Kel vin-Helmhol tz  billows  generated  in  the  laboratory  as  described 
by  Thorpe  (1973).  In  the  tank  experiments,  the  small-scale  turbulence  is 
generated  in  the  vortex  cores,  and  spreads  horizontally  to  amalgamate  the 
billows  at  dimensionless  time  of  roughly  x  =  x0  +  2,  where  x0  is  the  time  at 
which  the  fine-scale  turbulence  first  appears.  Since  t0  ~  3  in  the  numerical 
integration,  the  obsesrvation  that  the  turbulence  merges  horizontally  at  some 
time  between  x  =  4.3  and  x  =  5.8  is  quite  consistent  with  the  experimental 
data.  Furthermore,  the  thickness  of  the  region  of  turbulent  fluid  increases 
slowly  after  merging  in  the  horizontal  and  seems  to  reach  an  equilibrium  level 
near  the  end  of  the  integration.  This  is  again  consistent  with  the 

experimental  observations  of  Thorpe  (1973),  which  show  the  layer  increasing  in 
thickness  up  to  x  =  x0  +  7,  i.e.,  x  =  10  in  our  case.  Also  the 

non-dimensional  height,  Rj_  =  g(AT.h)/(2AU2)  where  h  is  the  thickness  of  the 
turbulent  layer,  reaches  a  value  of  about  0.4  in  accordance  with  the 


experiments. 


Some  insight  into  the  dynamics  is  obtainable  from  examination  of  the 
energy  budgets.  In  Figure  A. 3,  we  plot  the  roll  energy,  EK,  and  the 
small-scale  turbulence  energy,  EQ,  defined  by 


EK 
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}  dx  dz 


(A. 14) 


EQ  * 


(A. 16) 


where 


u  (z) 


(A. 16) 


We  can  see  the  rapid  initial  build-up  of  roll  energy  as  the  growing 
billow  extracts  kinetic  energy  from  the  mean  shear,  whilst  the  turbulence 
energy  EQ  remains  virtually  constant  at  its  initial  perturbation  level.  At 
roughly  \  -  2.5,  EQ  begins  to  increase  as  the  billow  starts  to  break,  and  EK 
levels  off  at  its  maximum  value.  The  kinetic  energy  is  effectively 
transferred  via  potential  energy  to  the  small-scale  turbulence  which  reaches  a 
peak  at  t  =  5.  The  large  scale  motion  loses  its  energy  very  rapidly  after 
breaking,  while  the  small-scale  turbulence  persists  throughout  the 
integration,  and  exhibits  a  relatively  slow  rate  of  decay.  Thorpe  (1973) 
shows  that  the  turbulence  persists  until  roughly  t  =  t0  +  12  =  15,  at  which 
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point  there  appears  to  be  a  collapse  of  the  turbulent  eddies  into  striated 
structures.  The  present  second-order  closure  model  of  the  small-scale  eddies 
provides  a  poor  representation  of  this  final  stage  of  the  decay  of  turbulence 
in  a  stably-stratified  medium,  and  we  have  therefore  not  attempted  to 
integrate  beyond  x  =  13. 

Profiles  of  the  initial  and  final  (x  =  11.8)  local  gradient  Richardson 
number , 


RiG 


g  9T/9Z 
T o  ( au/dz)^ 


(A. 17) 


are  shown  in  Figure  A. 4.  The  final  value  of  Rig  in  the  layer  is  not  quite 
constant,  but  lies  between  0.3  and  0.35.  This  is  again  consistent  with 
Thorpe's  estimate  of  the  final  Richardson  number.  This  value  is  significantly 
larger  than  the  critical  value  of  0.25,  so  that  the  mixing-process  apparently 
extracts  sufficient  energy  from  the  initially  unstable  flow  to  mix  the  mean 
profiles  further  than  the  point  at  which  no  more  energy  can  be  extracted.  The 
agreement  between  the  calculation  and  the  experiment  on  this  number  is  some 
confirmation  that  the  model  describes  the  energy  transfers  in  the  initial 
stages  reasonably  accurately.  The  basic  energy  source  is  kinetic  energy 
produced  by  the  billow  acting  on  the  mean  shear,  and  it  seems  likely  that  the 
final  amount  of  mixing  will  be  determined  by  the  amount  of  energy  that  the 
billow  can  extract  before  it  breaks.  If  the  turbulence  model  predicted  the 
wrong  time-scales  for  the  breaking  process,  we  would  be  unlikely  to  obtain  the 
correct  final  state. 

The  evolution  from  an  initial  Richardson  number  of  0.2  is  illustrated  in 
jures  A. 5  and  A. 6.  A  mesh  of  41  x  61  grid  points  was  used,  as  in  the 
previous  case,  with  a  domain  size  of  156  in  the  horizontal  again,  but  only  do 
in  the  vertical.  Contours  of  temperature  and  small-scale  turbulent  kinetic 
energy  at  x  =  3.0,  6.0,  9.0,  14.6,  and  20.6.  Clearly,  the  billow  is  much  more 
r  ic ted  in  the  vertical  due  to  the  fact  that  there  is  relatively  less 
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.4  -  Profiles  of  gradient  Richardson  number,  RiG,  for  the  Ri  =  0.1  case, 
ne  is  the  initial  profile  at  x  =  0,  dashed  line  is  the  profile  at  t  =  11 


T 


kinetic  energy  available  in  the  mean  shear  to  drive  the  instability  against 
the  restoring  buoyancy  forces.  The  higher  initial  Richardson  number  also 
makes  the  turbulent  kinetic  energy  production  in  the  braids  less  important,  so 
that  at  x  =  3  when  the  isotherms  are  only  just  being  overturned,  the  maximum 
energy  is  already  in  the  vortex  core.  Furthermore,  the  restricted  size  of  the 
vortex  core  appears  to  result  in  the  turbulence  being  spread  throughout  the 
region  rather  than  being  confined  to  the  outer  half  of  the  vortex  as  in  the 
Ri  =  0.1  case.  The  velocity  in  the  core  is  also  smaller  in  the  Ri  =  0.2  case, 
so  there  is  less  tendency  for  the  turbulence  to  be  swept  around  the  outer  part 
of  the  vortex. 

As  in  the  Ri  =  0.1  case,  the  turbulent  vortex  core  spreads  horizontally, 
with  small  scale  waves  running  along  the  top  and  bottom  of  the  layer  producing 
local  turbulence  maxima.  Finally,  the  mean  motion  is  totally  suppressed,  and 
the  isotherms  become  horizontal,  and  we  have  a  homogeneous  layer  of 
turbul ence. 

The  evolution  of  the  large-  and  small-scale  kinetic  energies  is  shown  in 
Figure  A. 7.  Comparing  this  with  Figure  A. 3,  we  see  that  the  large-scale 
roll-up  extracts  much  less  energy  than  the  Ri  =  0.1  case,  but  the  maximum  in 
small-scale  energy  is  almost  as  large  as  the  large-scale.  Since  the 
small-scale  energy  arises  from  the  potential  energy  created  by  the  roll-up 
process,  this  demonstrates  that  there  is  more  potential  energy  per  unit 
kinetic  energy  in  the  Ri  =  0.2  case.  This  should  not  be  a  surprise,  since  the 
Richardson  number  is  a  measure  of  the  ratio  of  potential  to  kinetic  energy  in 
the  initial  profile,  and  a  linearized  disturbance  will  contain  energy  in  the 
same  proportion. 

We  may  also  note  that  if  the  time  scale  for  the  Ri  =  0.2  case  is  divided 
by  two,  then  a  number  of  features  match  up  with  the  Ri  =0.1  case.  Defining 
x1  =  x/Ri ,  then  the  large  scale  energy  maximizes  at  x'  z  50  in  both  cases. 
Furthermore,  the  small-scale  energy  is  down  to  about  half  its  maximum  value  at 
x  =  100.  It  seems  that  the  time  scale  is  determined  by  the  shear  velocity  and 
the  layer  thickness,  while  the  buoyancy  time  scale  determines  the  details  of 
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the  roll-up  such  as  core  size  and  energy  partitioning. 

The  profile  of  gradient  Richardson  number  at  the  final  time,  t  =  20.6. 
is  shown  in  Figure  A. 8  for  the  case  with  initial  Ri  =  0.2.  The  final  state  is 
a  region  of  virtually  constant  Richardson  number  at  a  value  of  roughly  0.45, 
slightly  larger  than  the  case  with  Ri  =  0.1.  The  laboratory  experiments  of 
Thorpe  (1973)  indicate  that  the  final  Richardson  number  is  independent  of  the 
initial  value,  the  final  value  being  between  0.26  and  0.385,  but  this  result 
is  for  initial  Richardson  nimber  less  than  0.14. 

In  considering  atmospheric  observations,  we  must  recognize  the  fact  that 
the  high-resolution  radars  are  measuring  small-scale  refractive  index  changes. 
As  discussed  in  Chapter  II,  the  latter  are  very  closely  related  to  humidity 
and  temperature  fluctuations  rather  than  the  turbulence  energy  Hself.  Since 
humidity  is  a  scalar,  it  satisfies  the  same  equation  as  temperature,  and  it 
can  be  shown  that  the  mean  square  humidity  fluctuations  are  identical  to  the 
mean  square  temperature  fluctuation  provided  the  initial  profiles  are 
identical.  It  is  therefore  instructive  to  look  at  the  behavior  of  temperature 
fluctuations,  remembering  that  asymmetries  can  be  introduced  by  having  an 
asymmetric  humidity  profile.  In  fact,  the  length  scale  of  the 
energy-containing  turbulent  eddies  also  affects  the  refractive  index 
fluctuations  at  the  short  wavelengths  visible  to  the  radar,  but  variations  in 
the  length  scale  do  not  produce  any  significant  changes  in  the  pattern. 

Contours  of  the  mean  square  temperature  fluctuation,  la^MT2,  for  the  case 
with  initial  Ri  =  0.1  are  shown  in  Figure  A. 9  at  t  =  1.5,  2.8,  4.3,  and  8.8. 
The  most  obvious  feature  is  the  highly  localized  nature  of  "a7.  Initially,  "a7 
is  concentrated  in  the  developing  braids,  where  the  temperature  gradient  is 
very  large.  As  the  billow  begins  to  break,  there  are  local  patches  of  high 
temperature  variance  in  the  mixing  regions  around  the  edges  of  the  core,  but 
these  patches  do  not  spread  throughout  the  vortex  core  as  it  breaks.  Rather, 
as  the  vortex  breaks  down,  the  temperature  variance  vanishes,  and  only  remains 
in  the  braids  at  i  =  4.3.  As  the  braids  are  eroded  by  the  turbulent  mixing, 
the  temperature  variance  decays,  but  remains  in  narrow  bands  which  are  drawn 
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Isopleths  of  dimensionless  temperature  variance  62/aT2  for  the 
=  1.  (a)  t  =  1.5,  contour  interval  is  0.004. 


interval  is  0.02 


gure  A. 9(c)  -  t  =  4.3,  contour  interval  is  0.02. 


Figure  A. 9(d)  -  t  =  3.8,  contour  interval  is  0.002. 


out  horizontally  until  they  are  almost  flat.  This  qualitative  description 
also  applies  to  the  Ri  ■  0.2  case. 


The  difference  between  the  temperature  variance  and  velocity  variance  is 
also  evident  by  comparing  Figure  A. 10  with  Figures  A. 3  and  A. 7.  The 
integrated  value  of  the  temperature  variance  tends  to  both  grow  and  decay  more 
rapidly  then  the  integrated  velocity  variance. 

The  contrast  in  behavior  between  the  temperature  variance  and  turbulent 
kinetic  energy  may  be  somewhat  better  understood  by  examining  the  production 
terms  in  the  second-order  turbulence  correlations.  We  have 
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where  P(  )  denotes  production  terms.  The  only  other  terms  in  the  relevant 
equations  are  diffusion  and  dissipation  terms.  Now,  as  the  vortex  core  rolls 
up,  large  gradients  of  mean  quantities  are  generated,  and  consequently  the 
production  terms  are  also  large. 


Once  the  vortex  breaks  the  mean  gradients  are  rapidly  mixed  away,  and 
therefore  there  is  only  diffusion  and  dissipation  of  a7"  which  must  begin  to 
decay  immediately.  However,  there  is  still  production  of  vTe  due  to  the  term 
involving  the  temperature  variance.  Thus  heat  flux  does  not  decay  as  quickly 
as  "a7.  Furthermore,  wo  is  the  main  production  term  in  the  q2  equation  within 
the  core,  so  the  turbulent  kinetic  energy  continues  being  produced  while  a7  is 
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Figure  A. 10(a)  -  Evolution  of  the  total  volume-integrated  temperature  variance,  ET, 
for  the  case  with  Ri  =  0.1. 


0.45 


decaying.  Thus  q2  grows  in  the  core,  and  eventually  spreads  outward 
horizontally,  whilst  a7  vanishes  in  the  core,  and  only  remains  significant  in 
che  braids,  where  there  is  gradient  production.  When  the  braids  are  eroded, 
the  temperature  variance  is  slowly  diffused  and  dissipated,  while  the  mean 
shear  stretches  the  regions  of  high  "a7  into  nearly  horizontal  stratifications. 

The  differences  between  the  q2  and  "a7  patterns  are  thus  quite  plausible. 
We  note  that  this  dynamical  effect  could  not  be  reproduced  by  an  effective 
viscosity  model  which  has  no  equivalent  interaction  between  second-order 
quantities.  We  should  also  note  the  implication  that  high  resolution  radar 
returns  from  atmospheric  billow  events  of  this  type  should  always  show 
highly-structured  patterns  with  narrow  bands  of  high  intensity,  rather  than  a 
homogeneous  layer  of  turbulence.  This  is  quite  consistent  with  the 
observations  reported  in  the  literature  (See  Section  3 . 2  J .  FigureA.il 
showing  the  evolution  of  K-H  wave  breaking  as  received  by  radar  (Browning  and 
Watkins,  1970)  shows  the  same  features  as  exhibited  in  Figure  A. 9. 


A. 4  Sensitivity  Studies 

There  are  a  nunber  of  parameters  in  the  nunerical  integrations  presented 
in  the  previous  section  which  must  be  chosen  independently  from  the  physical 
constants  of  the  problem.  These  parameters  include  numerical  di screti zation 
quantities  such  as  the  grid  size  and  initial  conditions  for  the  dynamic 
variables.  We  do  not  include  the  empirical  constants  in  the  second-order 
turbulence  closure  model,  since  these  are  taken  to  be  fixed  from  comparison 
with  other  experimental  data.  It  is  important  to  obtain  some  indication  of 
the  dependence  of  the  results  or  these  external  parameters  if  we  are  to  gain  a 
useful  understanding  from  the  nunerical  integrations.  The  few  main  quantities 
which  need  to  be  chosen  are: 

A. 4.1  numerical  resolution 
A. 4. 2  initial  turbulence  level 

A. 4. 3  initial  turbulence  length  scale 

A. 4. 4  length  of  the  integration  domain 
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Figure  A. 11  -  Schematic  representation  of  the  life  cycle  of  an  individual 
Kelvin-Helmholtz  billow  based  on  the  data  in  the  earlier  figures.  Time 
progresses  from  right  to  left.  Thick  lines  correspond  to  the  detectable 
clear  air  radar  echo,  which  started  as  a  single  layer  at  1243  and  finished 
as  a  double  layer  at  1253  GMT.  Schematic  vertical  profiles  of  (A6//z)  are 
indicated  before  and  after  the  occurrence  of  Kelvin-Helmholtz  instability. 
(From  Browning  and  Watkins,  1970.) 


We  shall  investigate  variations  in  each  of  these  quantities  independently  for 
the  case  with  Ri  =  0.2. 

A. 4.1  Numerical  resolution.  The  integration  presented  in  the  previous 
section  was  carried  out  with  a  mesh  of  41  *  61  points;  this  was  chosen  to 
provide  adequate  resolution  of  the  small  scale  flow  features.  This  fact  is 
demonstrated  by  comparing  the  results  with  those  from  an  integration  using 
32  x  40  grid-points.  The  kinetic  energy  evolution  from  both  integrations  is 
shown  in  Figure  A. 12.  Unfortunately,  the  low  resolution  integration  had  a 
slightly  smaller  initial  vorticity  perturbation  which  is  responsible  for  some 
of  the  differences,  but  notwithstanding  this  initial  discrepancy,  the  overall 
budgets  are  within  about  10%. 

Contour  plots  from  the  two  runs  show  that  the  low  resolution  integration 
has  much  more  nunerical  noise  around  the  point  where  the  vortex  is  breaking 
and  gradients  are  highest.  However,  this  does  not  adversely  affect  the 
evolution,  as  can  be  seen  from  the  kinetic  energy  plots,  and  the  two 
integrations  do  remain  very  close  together  throughout  the  period. 

A. 4. 2  Initial  turbulence  level.  The  remaining  sensitivity  tests  are  all 
carried  out  using  the  low  resolution  model,  so  we  are  comparing  with  the  run 
in  Figure  A. 12(b) .  Figure  A. 13  shows  the  kinetic  energy  evolution  using  an 
initial  value  for  q^  of  3  x  10'^,  which  is  three  times  that  of  the  previous 
run.  The  main  effect  here  is  a  reduction  in  the  large  scale  energy  by  about 
12%,  and  an  increase  of  about  10%  in  the  small  scale  energy  maximum.  Thus  a 
higher  background  turbulence  level  retards  the  development  of  the  instability 
slightly,  and  also  allows  the  vortex  to  break  a  little  more  quickly, 
extracting  more  of  the  potential  energy  from  the  roll-up.  However,  the  effect 
is  not  large,  and  it  can  be  reasonably  assumed  that  the  results  are 
representative  of  small  background  turbulence  levels.  In  fact,  the  initial 
value  for  q2  was  chosen  as  large  as  reasonably  possible  in  order  to  provide 
some  turbulent  mixing  to  prevent  the  braids  from  becoming  too  thin  to  be 
resolved  by  the  numerical  grid.  Since  the  growth  of  q^  in  the  breaking  is  via 
the  buoyancy  term,  it  is  virtually  independent  of  the  initial  value  for  q2. 
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igure  A. 13  -  Evolution  of  eddy  kinetic  energies  for  Ri  =  0.2  with  d 
alues  of  q2/A(J2,  (a)  10"3. 


gure  A. 13(b) 


due  to  the  rapid  exponential  growth  in  the  initial  stages.  We  therefore  feel 
justified  in  choosing  q2  large  enough  to  control  the  braid  thickness,  as  long 
as  it  is  not  large  enough  to  influence  the  dynamics  of  the  breaking. 

Although  the  higher  value  of  initial  q2,  3  x  10_J,  appears  to  have  a 
small  effect  on  the  energy  budgets,  there  was  much  less  energy  in  the  short 
numerical  grid  modes,  so  this  run  was  chosen  as  the  base  integration  for 
comparison  of  turbulence  length  scale  and  domain  length  variations. 

A. 4. 3  Initial  turbulence  length  scale.  This  is  probably  the  most 
difficult  and  unfortunately  one  of  the  most  sensitive  parameters  to  choose. 
The  choice  of  initial  a  is  not  entirely  arbitrary,  since  it  must  bear  some 
relation  to  the  thickness  of  the  shear  layer,  6,  if  we  are  studying  turbulence 
generated  by  the  instabil ity  mechanism  itself.  We  chose  the  value  0.4o  for 
the  runs  presented  in  the  previous  section;  this  is  a  reasonable  value  for  a 
turbulence  length  scale,  but  cannot  be  justified  too  strongly.  In  order  to 
test  the  sensitivity,  we  carried  out  integrations  with  initial  A  =  0.2o  and 
0 . 8o ,  and  the  energy  evolutions  from  all  three  runs  are  shown  in  Figure  A. 14. 

Firstly,  the  large  scale  kinetic  energy  is  largely  unaffected  Dy  the 
choice  of  A.  The  only  difference  is  in  the  oscillation  in  the  decaying  phase 
cf  EK,  which  is  accentuated  by  small  values  of  a.  However,  the  small  scale 
energy  does  depend  quite  strongly  on  the  initial  value  of  a.  Changing  a  by  a 
factor  of  2  in  either  direction  from  0.4a  changes  the  maximum  turbulence 
energy  by  almost  a  factor  of  2  in  the  same  direction.  The  source  of  the 
problem  is  the  length  scale  equation,  which  does  not  allow  the  length  scale  to 
grow  sufficiently  quickly  during  the  breaking  stage,  so  that  the  turbulence 
does  remember  its  initial  condition  for  a  long  time.  The  length  scale  does 
grow  as  the  billow  breaks,  but  not  so  fast  as  q2  which  is  virtually 
independent  of  its  initial  condition. 

On  this  point,  all  we  are  able  to  say  is  that  0.46  is  a  reasonable 
initial  value  for  a,  and  the  results  of  this  section  illustrate  the  effect  of 
varying  a.  We  should  emphasize  that  the  qualitative  results,  both  in  the 
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energy  plots  and  the  details  of  the  contour  fields,  are  not  significantly 
changed  by  this  variation  in  A.  Apparently  the  initial  value  for  A  ought  to 
be  fixed  by  comparison  of  turbulence  levels  with  experimental  observation,  or 
an  improved  dynamical  equation  for  A  is  needed. 

A. 4. 4  Integration  domain  length.  We  choose  the  wavelength  of  the  unstable 
mode  we  are  studying  by  fixing  the  length  of  our  periodic  domain.  The  runs  in 
the  previous  section  used  a  length  of  156,  which  was  suggested  by  the  laminar 
numerical  integrations  of  Patnaik,  et  al.  (1976)  as  giving  the  mode  which  grew 
to  largest  amplitude. 

Figure  A. 15  shows  the  energies  from  three  runs  with  domain  lengths  of 
106,  156,  and  206.  Allowing  for  the  length  of  the  domain,  the  short  domain 
clearly  restricts  the  growth  of  the  billow  significantly;  the  large  scale 
energy  maximum  is  nearly  a  quarter  of  the  156  domain  case  (i.e.,  the  energy 
density  is  about  35%),  and  the  turbulence  energy  is  about  half  (i.e.,  energy 
density  is  70%).  It  is  therefore  unlikely  that  this  mode  would  be  the 
dominant  finite  amplitude  moae,  since  the  longer  modes  would  continue  to  grow. 
The  long  domain,  206,  has  a  significantly  larger  energy  in  the  vortex,  but 
allowing  for  the  length  of  the  domain,  the  energy  density  per  unit  horizontal 
length  is  only  a  little  longer.  The  turbulence  energy  per  unit  length  is 
actually  smaller. 

Thus  the  case  with  a  domain  length  of  156  seems  to  represent  the  most 
efficient  mode  for  producing  kinetic  energy,  although  there  is  not  a  great 
deal  of  difference  with  a  length  of  206. 

A. 5  Conclusions 

This  two-dimensional  numerical  study  of  breaking  Kel vin-Helmholtz  waves 
♦  using  a  second-order  turbulence  closure  model  to  describe  the  snail-scale 

turbulence  has  shown  agreement  with  experimental  observations  on  several 
points.  More  importantly,  given  the  lack  of  precise  quantitative  experimental 
data,  the  nunerical  integrations  have  suggested  several  general  features  of 
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.15  -  Evolution  of  eddy  kinetic  energies  for  Ri  =  0.2  with  di f ft 
enghts,  (a)  106. 
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the  dynamics  of  the  fine-scale  turbulence  in  addition  to  the  specific 
quantitative  predictions  of  the  particular  cases  studied.  These  general 
points  await  further  experimental  evidence  to  provide  confirmation  or 
refutation. 

The  results  demonstrate  the  utility  of  the  closure  model  in  studying  this 
problem;  closure  model  contains  most  of  the  important  dynamical  effects 
necessary  for  a  description  of  the  flow.  Furthermore,  it  is  probably  the 
minimum  sophistication  capable  of  describing  the  fine-scale  turbulence,  since 
it  has  shown  that  the  turbulence  develops  on  the  same  time  scale  as  the  mean 
flow,  so  that  equilibrium  assumptions  would  not  be  correct. 

The  numerical  results  indicate  the  correct  qualitative  behavior,  namely 
the  billow  breaks  to  produce  turbulence  in  the  vortex  core,  which  then  spreads 
horizontally  to  form  a  homogeneous  layer.  The  layer  also  spreads  a  little 
vertically.  The  time-scale  and  extents  of  the  spreading  are  in  broad 
agreement  with  the  visual  observations  of  Thorpe  (1973).  The  final  value  of 
the  Richardson  number  for  an  initial  Ri  of  0.1  is  consistent  with  Thorpe's 
measured  value  of  0.32  ±  0.06.  Moreover,  Thorpe  notes  no  dependence  on 
initial  Ri  for  Ri  <  0.14,  whilst  our  final  value  for  Ri  =  0.2  is  significantly 
higher.  It  remains  to  be  clarified  whether  the  experiments  would  show  such  a 
trend  at  higher  values  of  Ri. 

The  integrations  at  different  initial  Richardson  numbers  suggests  that 
the  advective  time  scale  j(,/aU,  where  i  is  the  wavelength  of  the  billow  and  au 
the  velociity  difference  across  the  shear  layer,  is  the  important  time  scale 
for  the  growth,  breaking,  and  subsequent  decay  of  the  turbulence.  The 
buoyancy  time  scale  affects  the  type  of  flow  since  the  ratio  of  the  time 
scales  is  proportional  to  the  square  root  of  the  Richardson  nunber.  It  is 
quite  possible  that  the  buoyancy  time  scale  is  also  important  in  the  late 
stages  of  the  decay,  when  the  turbulent  eddies  collapse  as  the  mean  density 
gradient  is  re-established.  We  have  not  attempted  to  seriously  study  this 
part  of  the  evolution,  since  the  closure  model  has  deficiencies  in  the  limit 
of  strongly-stratified  turbulence. 
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The  breaking  of  the  billow  itself  is  seen  to  be  generated  by  a  convective 
instability  of  the  overturned  fluid  at  the  top  and  bottom  of  the  vortex  core. 
This  turbulence  is  swept  around  the  vortex,  and  actually  remains  largely  in 
the  outer  part  of  the  vortex  for  a  significant  length  of  time  in  the  small 
Richardson  number  case.  The  turbulence  eventually  spreads  throughout  the 
layer,  and  amalgamates  to  form  a  horizontally  homogeneous  layer. 

Finally,  the  results  'ndicate  that  the  temperature  variance,  evolves 
quite  differently  from  the  turbulent-kinetic  energy.  Different  buoyancy 
production  terms  in  the  second-order  turbulence  correlation  equations  imply 
the  rapid  decay  of  in  the  breaking  vortex  core,  so  that  the  regions  of  high 
"a7  remain  highly  localized  in  space.  Effectively,  az  is  generated  initially 
in  the  braids  of  the  developing  billow,  and  these  long  thin  patches  are  then 
stretched  out  by  the  mean  shear  as  they  diffuse  and  decay  after  the  braids 
have  been  moved  away.  The  evolution  of  the  Kel vin-Helmhol tz  wave  breaking  as 
observed  by  radar  in  the  atmosphere  is  quite  consistent  with  the  evolution  of 
the  variance  of  humidity  or  temperature  predicted  by  the  model  calculation. 
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